From 444db1ae270daaec72e809f9077b6d4b95d0ecba Mon Sep 17 00:00:00 2001 From: Chris Geoga Date: Mon, 8 May 2023 10:39:35 -0400 Subject: [PATCH 01/17] Attempt at an extension/weakdep with an Enzyme rule. Not currently functional. --- LICENSE | 2 +- Project.toml | 8 +++++++- ext/EnzymeRules.jl | 50 ++++++++++++++++++++++++++++++++++++++++++++++ 3 files changed, 58 insertions(+), 2 deletions(-) create mode 100644 ext/EnzymeRules.jl diff --git a/LICENSE b/LICENSE index 51d9b9d..0bd745b 100644 --- a/LICENSE +++ b/LICENSE @@ -1,6 +1,6 @@ MIT License -Copyright (c) 2021-2022 Michael Helton, Oscar Smith, and contributors +Copyright (c) 2021-2022 Michael Helton, Oscar Smith, Chris Geoga, and contributors Permission is hereby granted, free of charge, to any person obtaining a copy of this software and associated documentation files (the "Software"), to deal diff --git a/Project.toml b/Project.toml index cb22803..b980beb 100644 --- a/Project.toml +++ b/Project.toml @@ -6,11 +6,17 @@ version = "0.3.0-DEV" SIMDMath = "5443be0b-e40a-4f70-a07e-dcd652efc383" [compat] -julia = "1.8" SIMDMath = "0.2.5" +julia = "1.8" [extras] Test = "8dfed614-e22c-5e08-85e1-65c5234f0b40" +[weakdeps] +Enzyme = "7da242da-08ed-463a-9acd-ee780be4f1d9" + +[extensions] +EnzymeRules = "Enzyme" + [targets] test = ["Test"] diff --git a/ext/EnzymeRules.jl b/ext/EnzymeRules.jl new file mode 100644 index 0000000..24b9a69 --- /dev/null +++ b/ext/EnzymeRules.jl @@ -0,0 +1,50 @@ + +# TODO (cg 2023/05/08 10:38): As of now, this throws the error +# +# ERROR; LoadError: UndefVarError: `forward` not defined +# +# which is weird, because this code works fine in its own scope and not as an +# extension. + +module EnzymeRules + + # TODO (cg 2023/05/08 10:02): Compat of any kind. + # + # TODO (cg 2023/05/08 10:27): This only works on the master branch of Enzyme. + # Which you can only use by doing + # ] add EnzymeCore#main + # ] add Enzyme#main + # so that's not great. + + using Bessels, Enzyme + using Enzyme.EnzymeRules + using Bessels.Math + + # A manual method that takes an NTuple of partial sum terms and checks if it is + # exactly converged before computing the Levin sequence transformation. If it is + # exactly converged, the sequence transformation will create a divide-by-zero + # problem. See + # + # https://github.com/JuliaMath/Bessels.jl/issues/96 + # + # and links with for discussion. + # + # TODO (cg 2023/05/08 10:00): I'm not entirely sure how best to "generalize" + # this to cases like a return type of DuplicatedNoNeed, or something being a + # `Enzyme.Const`. I have an open issue on the Enzyme.jl repo at + # + # https://github.com/EnzymeAD/Enzyme.jl/issues/786 + # + # that gets at this problem a bit. But it's a weird request and I'm sure Billy + # has a lot of asks on his time. + function EnzymeRules.forward(func::Const{typeof(levin_transform)}, + ::Type{<:Duplicated}, + s::Duplicated, + w::Duplicated) + (sv, dv, N) = (s.val, s.dval, length(s.val)) + ls = (sv[N-1] == sv[N]) ? sv[N] : levin_transform(sv, w.val) + dls = (dv[N-1] == dv[N]) ? dv[N] : levin_transform(dv, w.dval) + Duplicated(ls, dls) + end + +end From 17b90ca2d871ba5cf6c697ec56f8ad2ddb55bb59 Mon Sep 17 00:00:00 2001 From: Chris Geoga Date: Mon, 8 May 2023 10:54:10 -0400 Subject: [PATCH 02/17] Switch weakdep to EnzymeCore. --- Project.toml | 4 ++-- ext/EnzymeRules.jl | 10 ++++++---- 2 files changed, 8 insertions(+), 6 deletions(-) diff --git a/Project.toml b/Project.toml index b980beb..c83a4e3 100644 --- a/Project.toml +++ b/Project.toml @@ -13,10 +13,10 @@ julia = "1.8" Test = "8dfed614-e22c-5e08-85e1-65c5234f0b40" [weakdeps] -Enzyme = "7da242da-08ed-463a-9acd-ee780be4f1d9" +EnzymeCore = "f151be2c-9106-41f4-ab19-57ee4f262869" [extensions] -EnzymeRules = "Enzyme" +EnzymeRules = "EnzymeCore" [targets] test = ["Test"] diff --git a/ext/EnzymeRules.jl b/ext/EnzymeRules.jl index 24b9a69..c6bf56c 100644 --- a/ext/EnzymeRules.jl +++ b/ext/EnzymeRules.jl @@ -1,4 +1,3 @@ - # TODO (cg 2023/05/08 10:38): As of now, this throws the error # # ERROR; LoadError: UndefVarError: `forward` not defined @@ -16,8 +15,8 @@ module EnzymeRules # ] add Enzyme#main # so that's not great. - using Bessels, Enzyme - using Enzyme.EnzymeRules + using Bessels, EnzymeCore + using EnzymeCore.EnzymeRules using Bessels.Math # A manual method that takes an NTuple of partial sum terms and checks if it is @@ -31,7 +30,10 @@ module EnzymeRules # # TODO (cg 2023/05/08 10:00): I'm not entirely sure how best to "generalize" # this to cases like a return type of DuplicatedNoNeed, or something being a - # `Enzyme.Const`. I have an open issue on the Enzyme.jl repo at + # `Enzyme.Const`. These shouldn't in principle affect the "point" of this + # function (which is just to check for convergence before applying a + # function), but on its face this approach would mean I need a lot of + # hand-written extra methods. I have an open issue on the Enzyme.jl repo at # # https://github.com/EnzymeAD/Enzyme.jl/issues/786 # From 766cff5cf354041d7e40dcce7ababc505613f8e4 Mon Sep 17 00:00:00 2001 From: Chris Geoga Date: Mon, 8 May 2023 11:06:20 -0400 Subject: [PATCH 03/17] Add 2023 to license. --- LICENSE | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/LICENSE b/LICENSE index 0bd745b..feed8ff 100644 --- a/LICENSE +++ b/LICENSE @@ -1,6 +1,6 @@ MIT License -Copyright (c) 2021-2022 Michael Helton, Oscar Smith, Chris Geoga, and contributors +Copyright (c) 2021-2023 Michael Helton, Oscar Smith, Chris Geoga, and contributors Permission is hereby granted, free of charge, to any person obtaining a copy of this software and associated documentation files (the "Software"), to deal From 6db401efd1f53749a0a4e3256c599e50972be3b3 Mon Sep 17 00:00:00 2001 From: Chris Geoga Date: Mon, 8 May 2023 11:29:06 -0400 Subject: [PATCH 04/17] Add dummy test to trigger failure. --- test/Manifest.toml | 300 +++++++++++++++++++++++++++++--------------- test/Project.toml | 2 + test/enzyme_test.jl | 13 ++ 3 files changed, 211 insertions(+), 104 deletions(-) create mode 100644 test/enzyme_test.jl diff --git a/test/Manifest.toml b/test/Manifest.toml index 113c384..9ccd834 100644 --- a/test/Manifest.toml +++ b/test/Manifest.toml @@ -1,244 +1,336 @@ # This file is machine-generated - editing it directly is not advised -[[ArgTools]] +manifest_format = "2.0" + +[[deps.Adapt]] +deps = ["LinearAlgebra", "Requires"] +git-tree-sha1 = "cc37d689f599e8df4f464b2fa3870ff7db7492ef" +uuid = "79e6a3ab-5dfb-504d-930d-738a2a938a0e" +version = "3.6.1" + + [deps.Adapt.extensions] + AdaptStaticArraysExt = "StaticArrays" + + [deps.Adapt.weakdeps] + StaticArrays = "90137ffa-7385-5640-81b9-e52037218182" + +[[deps.ArgTools]] uuid = "0dad84c5-d112-42e6-8d28-ef12dabb789f" version = "1.1.1" -[[Artifacts]] +[[deps.Artifacts]] uuid = "56f22d72-fd6d-98f1-02f0-08ddc0907c33" -[[Base64]] +[[deps.Base64]] uuid = "2a0f44e3-6c83-55bd-87e4-b1978d98bd5f" -[[ChainRulesCore]] +[[deps.CEnum]] +git-tree-sha1 = "eb4cb44a499229b3b8426dcfb5dd85333951ff90" +uuid = "fa961155-64e5-5f13-b03f-caf6b980ea82" +version = "0.4.2" + +[[deps.ChainRulesCore]] deps = ["Compat", "LinearAlgebra", "SparseArrays"] -git-tree-sha1 = "4c26b4e9e91ca528ea212927326ece5918a04b47" +git-tree-sha1 = "c6d890a52d2c4d55d326439580c3b8d0875a77d9" uuid = "d360d2e6-b24c-11e9-a2a3-2a2ae2dbcce4" -version = "1.11.2" +version = "1.15.7" -[[ChangesOfVariables]] -deps = ["ChainRulesCore", "LinearAlgebra", "Test"] -git-tree-sha1 = "bf98fa45a0a4cee295de98d4c1462be26345b9a1" -uuid = "9e997f8a-9a97-42d5-a9f1-ce6bfc15e2c0" -version = "0.1.2" - -[[Compat]] -deps = ["Base64", "Dates", "DelimitedFiles", "Distributed", "InteractiveUtils", "LibGit2", "Libdl", "LinearAlgebra", "Markdown", "Mmap", "Pkg", "Printf", "REPL", "Random", "SHA", "Serialization", "SharedArrays", "Sockets", "SparseArrays", "Statistics", "Test", "UUIDs", "Unicode"] -git-tree-sha1 = "44c37b4636bc54afac5c574d2d02b625349d6582" +[[deps.Compat]] +deps = ["UUIDs"] +git-tree-sha1 = "7a60c856b9fa189eb34f5f8a6f6b5529b7942957" uuid = "34da2185-b29b-5c13-b0c7-acf172513d20" -version = "3.41.0" +version = "4.6.1" +weakdeps = ["Dates", "LinearAlgebra"] + + [deps.Compat.extensions] + CompatLinearAlgebraExt = "LinearAlgebra" -[[CompilerSupportLibraries_jll]] +[[deps.CompilerSupportLibraries_jll]] deps = ["Artifacts", "Libdl"] uuid = "e66e0078-7015-5450-92f7-15fbd957f2ae" -version = "0.5.2+0" +version = "1.0.2+0" -[[Dates]] +[[deps.Dates]] deps = ["Printf"] uuid = "ade2ca70-3891-5945-98fb-dc099432e06a" -[[DelimitedFiles]] -deps = ["Mmap"] -uuid = "8bb1440f-4735-579b-a4ab-409b98df4dab" - -[[Distributed]] -deps = ["Random", "Serialization", "Sockets"] -uuid = "8ba89e20-285c-5b6f-9357-94700520ee1b" - -[[DocStringExtensions]] +[[deps.DocStringExtensions]] deps = ["LibGit2"] -git-tree-sha1 = "b19534d1895d702889b219c382a6e18010797f0b" +git-tree-sha1 = "2fb1e02f2b635d0845df5d7c167fec4dd739b00d" uuid = "ffbed154-4ef7-542d-bbb7-c09d3a79fcae" -version = "0.8.6" +version = "0.9.3" -[[Downloads]] +[[deps.Downloads]] deps = ["ArgTools", "FileWatching", "LibCURL", "NetworkOptions"] uuid = "f43a241f-c20a-4ad4-852c-f6b1247861c6" version = "1.6.0" -[[FileWatching]] +[[deps.Enzyme]] +deps = ["CEnum", "EnzymeCore", "Enzyme_jll", "GPUCompiler", "LLVM", "Libdl", "LinearAlgebra", "ObjectFile", "Printf", "Random"] +git-tree-sha1 = "ad92a24f73bf8ffd1c23b2b5f59c9a761d7e9081" +repo-rev = "main" +repo-url = "https://github.com/EnzymeAD/Enzyme.jl.git" +uuid = "7da242da-08ed-463a-9acd-ee780be4f1d9" +version = "0.11.1" + +[[deps.EnzymeCore]] +deps = ["Adapt"] +git-tree-sha1 = "d0840cfff51e34729d20fd7d0a13938dc983878b" +repo-rev = "main" +repo-subdir = "lib/EnzymeCore" +repo-url = "https://github.com/EnzymeAD/Enzyme.jl.git" +uuid = "f151be2c-9106-41f4-ab19-57ee4f262869" +version = "0.3.0" + +[[deps.Enzyme_jll]] +deps = ["Artifacts", "JLLWrappers", "LazyArtifacts", "Libdl", "TOML"] +git-tree-sha1 = "5aaf389552c3560da95e2767080a22037eed90a5" +uuid = "7cc45869-7501-5eee-bdea-0790c847d4ef" +version = "0.0.59+0" + +[[deps.ExprTools]] +git-tree-sha1 = "c1d06d129da9f55715c6c212866f5b1bddc5fa00" +uuid = "e2ba6199-217a-4e67-a87a-7c52f15ade04" +version = "0.1.9" + +[[deps.FileWatching]] uuid = "7b1f6079-737a-58dc-b8bc-7a2ca5c1b5ee" -[[InteractiveUtils]] +[[deps.GPUCompiler]] +deps = ["ExprTools", "InteractiveUtils", "LLVM", "Libdl", "Logging", "Scratch", "TimerOutputs", "UUIDs"] +git-tree-sha1 = "e9a9173cd77e16509cdf9c1663fda19b22a518b7" +uuid = "61eb1bfa-7361-4325-ad38-22787b887f55" +version = "0.19.3" + +[[deps.InteractiveUtils]] deps = ["Markdown"] uuid = "b77e0a4c-d291-57a0-90e8-8db25a27a240" -[[InverseFunctions]] -deps = ["Test"] -git-tree-sha1 = "a7254c0acd8e62f1ac75ad24d5db43f5f19f3c65" -uuid = "3587e190-3f89-42d0-90ee-14403ec27112" -version = "0.1.2" - -[[IrrationalConstants]] +[[deps.IrrationalConstants]] git-tree-sha1 = "7fd44fd4ff43fc60815f8e764c0f352b83c49151" uuid = "92d709cd-6900-40b7-9082-c6be49f344b6" version = "0.1.1" -[[JLLWrappers]] +[[deps.JLLWrappers]] deps = ["Preferences"] -git-tree-sha1 = "642a199af8b68253517b80bd3bfd17eb4e84df6e" +git-tree-sha1 = "abc9885a7ca2052a736a600f7fa66209f96506e1" uuid = "692b3bcd-3c85-4b1f-b108-f13ce0eb3210" -version = "1.3.0" +version = "1.4.1" -[[LibCURL]] +[[deps.LLVM]] +deps = ["CEnum", "LLVMExtra_jll", "Libdl", "Printf", "Unicode"] +git-tree-sha1 = "a8960cae30b42b66dd41808beb76490519f6f9e2" +uuid = "929cbde3-209d-540e-8aea-75f648917ca0" +version = "5.0.0" + +[[deps.LLVMExtra_jll]] +deps = ["Artifacts", "JLLWrappers", "LazyArtifacts", "Libdl", "TOML"] +git-tree-sha1 = "09b7505cc0b1cee87e5d4a26eea61d2e1b0dcd35" +uuid = "dad2f222-ce93-54a1-a47d-0025e8a3acab" +version = "0.0.21+0" + +[[deps.LazyArtifacts]] +deps = ["Artifacts", "Pkg"] +uuid = "4af54fe1-eca0-43a8-85a7-787d91b784e3" + +[[deps.LibCURL]] deps = ["LibCURL_jll", "MozillaCACerts_jll"] uuid = "b27032c2-a3e7-50c8-80cd-2d36dbcbfd21" version = "0.6.3" -[[LibCURL_jll]] +[[deps.LibCURL_jll]] deps = ["Artifacts", "LibSSH2_jll", "Libdl", "MbedTLS_jll", "Zlib_jll", "nghttp2_jll"] uuid = "deac9b47-8bc7-5906-a0fe-35ac56dc84c0" version = "7.84.0+0" -[[LibGit2]] +[[deps.LibGit2]] deps = ["Base64", "NetworkOptions", "Printf", "SHA"] uuid = "76f85450-5226-5b5a-8eaa-529ad045b433" -[[LibSSH2_jll]] +[[deps.LibSSH2_jll]] deps = ["Artifacts", "Libdl", "MbedTLS_jll"] uuid = "29816b5a-b9ab-546f-933c-edad1886dfa8" version = "1.10.2+0" -[[Libdl]] +[[deps.Libdl]] uuid = "8f399da3-3557-5675-b5ff-fb832c97cbdb" -[[LinearAlgebra]] -deps = ["Libdl", "libblastrampoline_jll"] +[[deps.LinearAlgebra]] +deps = ["Libdl", "OpenBLAS_jll", "libblastrampoline_jll"] uuid = "37e2e46d-f89d-539d-b4ee-838fcccc9c8e" -[[LogExpFunctions]] -deps = ["ChainRulesCore", "ChangesOfVariables", "DocStringExtensions", "InverseFunctions", "IrrationalConstants", "LinearAlgebra"] -git-tree-sha1 = "be9eef9f9d78cecb6f262f3c10da151a6c5ab827" +[[deps.LogExpFunctions]] +deps = ["DocStringExtensions", "IrrationalConstants", "LinearAlgebra"] +git-tree-sha1 = "0a1b7c2863e44523180fdb3146534e265a91870b" uuid = "2ab3a3ac-af41-5b50-aa03-7779005ae688" -version = "0.3.5" +version = "0.3.23" + + [deps.LogExpFunctions.extensions] + LogExpFunctionsChainRulesCoreExt = "ChainRulesCore" + LogExpFunctionsChangesOfVariablesExt = "ChangesOfVariables" + LogExpFunctionsInverseFunctionsExt = "InverseFunctions" -[[Logging]] + [deps.LogExpFunctions.weakdeps] + ChainRulesCore = "d360d2e6-b24c-11e9-a2a3-2a2ae2dbcce4" + ChangesOfVariables = "9e997f8a-9a97-42d5-a9f1-ce6bfc15e2c0" + InverseFunctions = "3587e190-3f89-42d0-90ee-14403ec27112" + +[[deps.Logging]] uuid = "56ddb016-857b-54e1-b83d-db4d58db5568" -[[Markdown]] +[[deps.Markdown]] deps = ["Base64"] uuid = "d6f4376e-aef5-505a-96c1-9c027394607a" -[[MbedTLS_jll]] +[[deps.MbedTLS_jll]] deps = ["Artifacts", "Libdl"] uuid = "c8ffd9c3-330d-5841-b78e-0817d7145fa1" -version = "2.28.0+0" - -[[Mmap]] -uuid = "a63ad114-7e13-5084-954f-fe012c677804" +version = "2.28.2+0" -[[MozillaCACerts_jll]] +[[deps.MozillaCACerts_jll]] uuid = "14a3606d-f60d-562e-9121-12d972cd8159" -version = "2022.2.1" +version = "2022.10.11" -[[NetworkOptions]] +[[deps.NetworkOptions]] uuid = "ca575930-c2e3-43a9-ace4-1e988b2c1908" version = "1.2.0" -[[OpenBLAS_jll]] +[[deps.ObjectFile]] +deps = ["Reexport", "StructIO"] +git-tree-sha1 = "55ce61d43409b1fb0279d1781bf3b0f22c83ab3b" +uuid = "d8793406-e978-5875-9003-1fc021f44a92" +version = "0.3.7" + +[[deps.OpenBLAS_jll]] deps = ["Artifacts", "CompilerSupportLibraries_jll", "Libdl"] uuid = "4536629a-c528-5b80-bd46-f80d51c5b363" -version = "0.3.20+0" +version = "0.3.21+4" -[[OpenLibm_jll]] +[[deps.OpenLibm_jll]] deps = ["Artifacts", "Libdl"] uuid = "05823500-19ac-5b8b-9628-191a04bc5112" version = "0.8.1+0" -[[OpenSpecFun_jll]] +[[deps.OpenSpecFun_jll]] deps = ["Artifacts", "CompilerSupportLibraries_jll", "JLLWrappers", "Libdl", "Pkg"] git-tree-sha1 = "13652491f6856acfd2db29360e1bbcd4565d04f1" uuid = "efe28fd5-8261-553b-a9e1-b2916fc3738e" version = "0.5.5+0" -[[Pkg]] -deps = ["Artifacts", "Dates", "Downloads", "LibGit2", "Libdl", "Logging", "Markdown", "Printf", "REPL", "Random", "SHA", "Serialization", "TOML", "Tar", "UUIDs", "p7zip_jll"] +[[deps.Pkg]] +deps = ["Artifacts", "Dates", "Downloads", "FileWatching", "LibGit2", "Libdl", "Logging", "Markdown", "Printf", "REPL", "Random", "SHA", "Serialization", "TOML", "Tar", "UUIDs", "p7zip_jll"] uuid = "44cfe95a-1eb2-52ea-b672-e2afdf69b78f" -version = "1.8.0" +version = "1.9.0" -[[Preferences]] +[[deps.Preferences]] deps = ["TOML"] -git-tree-sha1 = "00cfd92944ca9c760982747e9a1d0d5d86ab1e5a" +git-tree-sha1 = "7eb1686b4f04b82f96ed7a4ea5890a4f0c7a09f1" uuid = "21216c6a-2e73-6563-6e65-726566657250" -version = "1.2.2" +version = "1.4.0" -[[Printf]] +[[deps.Printf]] deps = ["Unicode"] uuid = "de0858da-6303-5e67-8744-51eddeeeb8d7" -[[REPL]] +[[deps.REPL]] deps = ["InteractiveUtils", "Markdown", "Sockets", "Unicode"] uuid = "3fa0cd96-eef1-5676-8a61-b3b8758bbffb" -[[Random]] +[[deps.Random]] deps = ["SHA", "Serialization"] uuid = "9a3f8284-a2c9-5f02-9a11-845980a1fd5c" -[[SHA]] +[[deps.Reexport]] +git-tree-sha1 = "45e428421666073eab6f2da5c9d310d99bb12f9b" +uuid = "189a3867-3050-52da-a836-e630ba90ab69" +version = "1.2.2" + +[[deps.Requires]] +deps = ["UUIDs"] +git-tree-sha1 = "838a3a4188e2ded87a4f9f184b4b0d78a1e91cb7" +uuid = "ae029012-a4dd-5104-9daa-d747884805df" +version = "1.3.0" + +[[deps.SHA]] uuid = "ea8e919c-243c-51af-8825-aaa63cd721ce" version = "0.7.0" -[[Serialization]] -uuid = "9e88b42a-f829-5b0c-bbe9-9e923198166b" +[[deps.Scratch]] +deps = ["Dates"] +git-tree-sha1 = "30449ee12237627992a99d5e30ae63e4d78cd24a" +uuid = "6c6a2e73-6563-6170-7368-637461726353" +version = "1.2.0" -[[SharedArrays]] -deps = ["Distributed", "Mmap", "Random", "Serialization"] -uuid = "1a1011a3-84de-559e-8e89-a11a2f7dc383" +[[deps.Serialization]] +uuid = "9e88b42a-f829-5b0c-bbe9-9e923198166b" -[[Sockets]] +[[deps.Sockets]] uuid = "6462fe0b-24de-5631-8697-dd941f90decc" -[[SparseArrays]] -deps = ["LinearAlgebra", "Random"] +[[deps.SparseArrays]] +deps = ["Libdl", "LinearAlgebra", "Random", "Serialization", "SuiteSparse_jll"] uuid = "2f01184e-e22b-5df5-ae63-d93ebab69eaf" -[[SpecialFunctions]] +[[deps.SpecialFunctions]] deps = ["ChainRulesCore", "IrrationalConstants", "LogExpFunctions", "OpenLibm_jll", "OpenSpecFun_jll"] git-tree-sha1 = "e08890d19787ec25029113e88c34ec20cac1c91e" uuid = "276daf66-3868-5448-9aa4-cd146d93841b" version = "2.0.0" -[[Statistics]] -deps = ["LinearAlgebra", "SparseArrays"] -uuid = "10745b16-79ce-11e8-11f9-7d13ad32a3b2" +[[deps.StructIO]] +deps = ["Test"] +git-tree-sha1 = "010dc73c7146869c042b49adcdb6bf528c12e859" +uuid = "53d494c1-5632-5724-8f4c-31dff12d585f" +version = "0.3.0" + +[[deps.SuiteSparse_jll]] +deps = ["Artifacts", "Libdl", "Pkg", "libblastrampoline_jll"] +uuid = "bea87d4a-7f5b-5778-9afe-8cc45184846c" +version = "5.10.1+6" -[[TOML]] +[[deps.TOML]] deps = ["Dates"] uuid = "fa267f1f-6049-4f14-aa54-33bafae1ed76" -version = "1.0.0" +version = "1.0.3" -[[Tar]] +[[deps.Tar]] deps = ["ArgTools", "SHA"] uuid = "a4e569a6-e804-4fa4-b0f3-eef7a1d5b13e" -version = "1.10.1" +version = "1.10.0" -[[Test]] +[[deps.Test]] deps = ["InteractiveUtils", "Logging", "Random", "Serialization"] uuid = "8dfed614-e22c-5e08-85e1-65c5234f0b40" -[[UUIDs]] +[[deps.TimerOutputs]] +deps = ["ExprTools", "Printf"] +git-tree-sha1 = "f548a9e9c490030e545f72074a41edfd0e5bcdd7" +uuid = "a759f4b9-e2f1-59dc-863e-4aeb61b1ea8f" +version = "0.5.23" + +[[deps.UUIDs]] deps = ["Random", "SHA"] uuid = "cf7118a7-6976-5b1a-9a39-7adc72f591a4" -[[Unicode]] +[[deps.Unicode]] uuid = "4ec0a83e-493e-50e2-b9ac-8f72acf5a8f5" -[[Zlib_jll]] +[[deps.Zlib_jll]] deps = ["Libdl"] uuid = "83775a58-1f1d-513f-b197-d71354ab007a" -version = "1.2.12+3" +version = "1.2.13+0" -[[libblastrampoline_jll]] -deps = ["Artifacts", "Libdl", "OpenBLAS_jll"] +[[deps.libblastrampoline_jll]] +deps = ["Artifacts", "Libdl"] uuid = "8e850b90-86db-534c-a0d3-1478176c7d93" -version = "5.1.1+0" +version = "5.7.0+0" -[[nghttp2_jll]] +[[deps.nghttp2_jll]] deps = ["Artifacts", "Libdl"] uuid = "8e850ede-7688-5339-a07c-302acd2aaf8d" version = "1.48.0+0" -[[p7zip_jll]] +[[deps.p7zip_jll]] deps = ["Artifacts", "Libdl"] uuid = "3f19e933-33d8-53b3-aaab-bd5110c3b7a0" version = "17.4.0+0" diff --git a/test/Project.toml b/test/Project.toml index e7afe0b..b19b398 100644 --- a/test/Project.toml +++ b/test/Project.toml @@ -1,3 +1,5 @@ [deps] +Enzyme = "7da242da-08ed-463a-9acd-ee780be4f1d9" +EnzymeCore = "f151be2c-9106-41f4-ab19-57ee4f262869" SpecialFunctions = "276daf66-3868-5448-9aa4-cd146d93841b" Test = "8dfed614-e22c-5e08-85e1-65c5234f0b40" diff --git a/test/enzyme_test.jl b/test/enzyme_test.jl new file mode 100644 index 0000000..533308c --- /dev/null +++ b/test/enzyme_test.jl @@ -0,0 +1,13 @@ + +# This is just a placeholder at the moment to try and make CI trigger the error. +# This should hit the special method created in the EnzymeRules extension/weakdep. + +using EnzymeCore, Enzyme +using Bessels.BesselFunctions +using BesselFunctions: besselkx_levin + +dbesselk_dv(v, x) = autodiff(Forward, _v->besselkx_levin(_v, x, Val(20)), + Duplicated, Duplicated(v, 1.0)) + +@assert isapprox(dbesselk_dv(1.5, 13.1), 0.0410929, atol=1e-5) + From dbb3821e89137e34a850cb1f63e9e1cf48e62519 Mon Sep 17 00:00:00 2001 From: Chris Geoga Date: Mon, 8 May 2023 11:39:12 -0400 Subject: [PATCH 05/17] Switching to tagged releases of Enzyme and EnzymeCore to fix the mysterious "cannot merge projects" error. --- test/Manifest.toml | 13 +++++-------- test/runtests.jl | 1 + 2 files changed, 6 insertions(+), 8 deletions(-) diff --git a/test/Manifest.toml b/test/Manifest.toml index 9ccd834..1667257 100644 --- a/test/Manifest.toml +++ b/test/Manifest.toml @@ -1,6 +1,8 @@ # This file is machine-generated - editing it directly is not advised +julia_version = "1.9.0-rc3" manifest_format = "2.0" +project_hash = "d08b9c18a9b1978d08c223ae5573bfa705b3affe" [[deps.Adapt]] deps = ["LinearAlgebra", "Requires"] @@ -67,26 +69,21 @@ version = "1.6.0" [[deps.Enzyme]] deps = ["CEnum", "EnzymeCore", "Enzyme_jll", "GPUCompiler", "LLVM", "Libdl", "LinearAlgebra", "ObjectFile", "Printf", "Random"] -git-tree-sha1 = "ad92a24f73bf8ffd1c23b2b5f59c9a761d7e9081" -repo-rev = "main" -repo-url = "https://github.com/EnzymeAD/Enzyme.jl.git" +git-tree-sha1 = "4478f8bf24785d9eabe09044549b0e81b9e12d68" uuid = "7da242da-08ed-463a-9acd-ee780be4f1d9" version = "0.11.1" [[deps.EnzymeCore]] deps = ["Adapt"] git-tree-sha1 = "d0840cfff51e34729d20fd7d0a13938dc983878b" -repo-rev = "main" -repo-subdir = "lib/EnzymeCore" -repo-url = "https://github.com/EnzymeAD/Enzyme.jl.git" uuid = "f151be2c-9106-41f4-ab19-57ee4f262869" version = "0.3.0" [[deps.Enzyme_jll]] deps = ["Artifacts", "JLLWrappers", "LazyArtifacts", "Libdl", "TOML"] -git-tree-sha1 = "5aaf389552c3560da95e2767080a22037eed90a5" +git-tree-sha1 = "b0f72433c4679db4df05c999f200d60cb78d1a27" uuid = "7cc45869-7501-5eee-bdea-0790c847d4ef" -version = "0.0.59+0" +version = "0.0.57+0" [[deps.ExprTools]] git-tree-sha1 = "c1d06d129da9f55715c6c212866f5b1bddc5fa00" diff --git a/test/runtests.jl b/test/runtests.jl index 69e6085..db22996 100644 --- a/test/runtests.jl +++ b/test/runtests.jl @@ -10,3 +10,4 @@ import SpecialFunctions @time @testset "gamma" begin include("gamma_test.jl") end @time @testset "airy" begin include("airy_test.jl") end @time @testset "sphericalbessel" begin include("sphericalbessel_test.jl") end +@time @testset "enzyme autodiff" begin include("enzyme_test.jl") end From e2d1f41441fe31007e4b9710d478d1636c5d9fe0 Mon Sep 17 00:00:00 2001 From: Chris Geoga Date: Tue, 9 May 2023 10:00:47 -0400 Subject: [PATCH 06/17] Extension fixed, test passing. --- Project.toml | 2 +- ext/{EnzymeRules.jl => BesselsEnzymeCoreExt.jl} | 15 +-------------- src/BesselFunctions/besselk.jl | 15 ++++++++------- src/Math/Math.jl | 5 +++-- test/enzyme_test.jl | 7 +++---- 5 files changed, 16 insertions(+), 28 deletions(-) rename ext/{EnzymeRules.jl => BesselsEnzymeCoreExt.jl} (78%) diff --git a/Project.toml b/Project.toml index c83a4e3..edd7ab9 100644 --- a/Project.toml +++ b/Project.toml @@ -16,7 +16,7 @@ Test = "8dfed614-e22c-5e08-85e1-65c5234f0b40" EnzymeCore = "f151be2c-9106-41f4-ab19-57ee4f262869" [extensions] -EnzymeRules = "EnzymeCore" +BesselsEnzymeCoreExt = "EnzymeCore" [targets] test = ["Test"] diff --git a/ext/EnzymeRules.jl b/ext/BesselsEnzymeCoreExt.jl similarity index 78% rename from ext/EnzymeRules.jl rename to ext/BesselsEnzymeCoreExt.jl index c6bf56c..d21493c 100644 --- a/ext/EnzymeRules.jl +++ b/ext/BesselsEnzymeCoreExt.jl @@ -1,19 +1,6 @@ -# TODO (cg 2023/05/08 10:38): As of now, this throws the error -# -# ERROR; LoadError: UndefVarError: `forward` not defined -# -# which is weird, because this code works fine in its own scope and not as an -# extension. - -module EnzymeRules +module BesselsEnzymeCoreExt # TODO (cg 2023/05/08 10:02): Compat of any kind. - # - # TODO (cg 2023/05/08 10:27): This only works on the master branch of Enzyme. - # Which you can only use by doing - # ] add EnzymeCore#main - # ] add Enzyme#main - # so that's not great. using Bessels, EnzymeCore using EnzymeCore.EnzymeRules diff --git a/src/BesselFunctions/besselk.jl b/src/BesselFunctions/besselk.jl index 625aab9..7407464 100644 --- a/src/BesselFunctions/besselk.jl +++ b/src/BesselFunctions/besselk.jl @@ -578,15 +578,16 @@ end @generated function besselkx_levin(v, x::T, ::Val{N}) where {T <: FloatTypes, N} :( begin - s_0 = zero(T) + s = zero(T) t = one(T) @nexprs $N i -> begin - s_{i} = s_{i-1} + t - t *= (4*v^2 - (2i - 1)^2) / (8 * x * i) - w_{i} = 1 / t - end - sequence = @ntuple $N i -> s_{i} - weights = @ntuple $N i -> w_{i} + s += t + t *= (4*v^2 - (2i - 1)^2) / (8 * x * i) + s_{i} = s + w_{i} = t + end + sequence = @ntuple $N i -> s_{i} + weights = @ntuple $N i -> w_{i} return levin_transform(sequence, weights) * sqrt(π / 2x) end ) diff --git a/src/Math/Math.jl b/src/Math/Math.jl index 8b1c840..386ca64 100644 --- a/src/Math/Math.jl +++ b/src/Math/Math.jl @@ -131,11 +131,12 @@ end #@inline levin_scale(B::T, n, k) where T = -(B + n) * (B + n + k)^(k - one(T)) / (B + n + k + one(T))^k @inline levin_scale(B::T, n, k) where T = -(B + n + k) * (B + n + k - 1) / ((B + n + 2k) * (B + n + 2k - 1)) -@inline @generated function levin_transform(s::NTuple{N, T}, w::NTuple{N, T}) where {N, T <: FloatTypes} +@inline @generated function levin_transform(s::NTuple{N, T}, + w::NTuple{N, T}) where {N, T <: FloatTypes} len = N - 1 :( begin - @nexprs $N i -> a_{i} = Vec{2, T}((s[i] * w[i], w[i])) + @nexprs $N i -> a_{i} = Vec{2, T}((s[i] / w[i], 1 / w[i])) @nexprs $len k -> (@nexprs ($len-k) i -> a_{i} = fmadd(a_{i}, levin_scale(one(T), i, k-1), a_{i+1})) return (a_1[1] / a_1[2]) end diff --git a/test/enzyme_test.jl b/test/enzyme_test.jl index 533308c..fcc2675 100644 --- a/test/enzyme_test.jl +++ b/test/enzyme_test.jl @@ -3,11 +3,10 @@ # This should hit the special method created in the EnzymeRules extension/weakdep. using EnzymeCore, Enzyme -using Bessels.BesselFunctions -using BesselFunctions: besselkx_levin +import Bessels.BesselFunctions: besselkx_levin dbesselk_dv(v, x) = autodiff(Forward, _v->besselkx_levin(_v, x, Val(20)), - Duplicated, Duplicated(v, 1.0)) + Duplicated, Duplicated(v, 1.0))[2] -@assert isapprox(dbesselk_dv(1.5, 13.1), 0.0410929, atol=1e-5) +@test isapprox(dbesselk_dv(1.5, 13.1), 0.0410929, atol=1e-5) From fc1e4a7fa594cd1fe409cd00320fbced966111a1 Mon Sep 17 00:00:00 2001 From: Chris Geoga Date: Tue, 9 May 2023 11:35:32 -0400 Subject: [PATCH 07/17] Test framework in place, but I don't trust these Arb derivatives. --- test/data/besselk/enzyme/README.md | 31 ++ .../enzyme/besselkx_levin_enzyme_tests.csv | 480 ++++++++++++++++++ test/enzyme_test.jl | 18 +- test/runtests.jl | 16 +- 4 files changed, 530 insertions(+), 15 deletions(-) create mode 100644 test/data/besselk/enzyme/README.md create mode 100644 test/data/besselk/enzyme/besselkx_levin_enzyme_tests.csv diff --git a/test/data/besselk/enzyme/README.md b/test/data/besselk/enzyme/README.md new file mode 100644 index 0000000..33bcb9f --- /dev/null +++ b/test/data/besselk/enzyme/README.md @@ -0,0 +1,31 @@ +```julia +using FiniteDifferences, ArbNumerics, DelimitedFiles + +if !(@isdefined vx) + const FD = central_fdm(10,1) + vgrid = sort(vcat(range(0.05, 15.0, length=20), [0.5, 1.5, 2.5, 3.5])) + xgrid = range(15.0, 30.0, length=20) + const vx = vec(collect(Iterators.product(vgrid, xgrid))) +end + +ArbNumerics.setprecision(ArbReal, digits=100) + +function arb_besselkx(v,x) + (av, ax) = ArbReal.((v,x)) + ArbNumerics.besselk(av, ax)*ArbNumerics.exp(ax) +end + +ref_values = map(vx) do vxj + (v,x) = vxj + dx = FD(_x->arb_besselkx(v, _x), x) + dv = FD(_v->arb_besselkx(_v, x), v) + Float64.((dv, dx)) +end + +out_matrix = hcat(getindex.(vx, 1), # test v argument + getindex.(vx, 2), # test x argument + getindex.(ref_values, 1), # test d/dv value + getindex.(ref_values, 2)) # test d/dx value + +writedlm("besselkx_levin_enzyme_tests.csv", out_matrix) +``` diff --git a/test/data/besselk/enzyme/besselkx_levin_enzyme_tests.csv b/test/data/besselk/enzyme/besselkx_levin_enzyme_tests.csv new file mode 100644 index 0000000..4cdc411 --- /dev/null +++ b/test/data/besselk/enzyme/besselkx_levin_enzyme_tests.csv @@ -0,0 +1,480 @@ +0.05 15.0 0.0010367551275860735 -0.010535444408076956 +0.5 15.0 0.01044909892351376 -0.010786835957018058 +0.8368421052631579 15.0 0.017739559706085498 -0.011253628340044715 +1.5 15.0 0.03339462144739582 -0.01294415009942762 +1.6236842105263158 15.0 0.03659268674207609 -0.013387189370521184 +2.4105263157894736 15.0 0.060052964058098846 -0.017373888634610528 +2.5 15.0 0.06315118081428582 -0.017978120147579314 +3.1973684210526314 15.0 0.09153794927366792 -0.024071625383024185 +3.5 15.0 0.10679821808137248 -0.027662120748450387 +3.9842105263157896 15.0 0.1362234972398303 -0.03500862169499547 +4.771052631578947 15.0 0.20238301287210217 -0.05287136846744076 +5.557894736842106 15.0 0.30374695994764694 -0.08239500819064587 +6.344736842105263 15.0 0.4636380330289938 -0.13208397871236013 +7.131578947368421 15.0 0.7228219860633713 -0.21742481203781536 +7.918421052631579 15.0 1.1539993886345965 -0.36729242248277877 +8.705263157894738 15.0 1.8899224810295472 -0.6365147679147345 +9.492105263157894 15.0 3.178648748030142 -1.1317368209637213 +10.278947368421052 15.0 5.4925111585656206 -2.064343411173481 +11.06578947368421 15.0 9.755179997903623 -3.862762417806891 +11.852631578947369 15.0 17.809093701857034 -7.414088150681344 +12.639473684210527 15.0 33.42266223980047 -14.596901632387276 +13.426315789473684 15.0 64.44838838752855 -29.471567616652393 +14.213157894736842 15.0 127.68673796168575 -61.01851035096151 +15.0 15.0 259.78299436259425 -129.48144922416114 +0.05 15.789473684210526 0.0009618017615888898 -0.009765556568247924 +0.5 15.789473684210526 0.009690053231145072 -0.009988315614370639 +0.8368421052631579 15.789473684210526 0.016439733325060442 -0.010399134668978188 +1.5 15.789473684210526 0.03087566747671769 -0.011886188577293443 +1.6236842105263158 15.789473684210526 0.03381256980329193 -0.0122738720493344 +2.4105263157894736 15.789473684210526 0.05522747048187236 -0.015757523702255714 +2.5 15.789473684210526 0.058038666965933475 -0.016282009216853487 +3.1973684210526314 15.789473684210526 0.08362918041182572 -0.02155695614180951 +3.5 15.789473684210526 0.09727720780815732 -0.024645549572446225 +3.9842105263157896 15.789473684210526 0.12342043461990765 -0.03092882160405939 +4.771052631578947 15.789473684210526 0.18152361114544602 -0.046042595633940836 +5.557894736842106 15.789473684210526 0.2692415361359061 -0.07068159896315275 +6.344736842105263 15.789473684210526 0.4055175223058745 -0.11150082248932047 +7.131578947368421 15.789473684210526 0.622835091385322 -0.18046448445800095 +7.918421052631579 15.789473684210526 0.9782051119542616 -0.29944048626010167 +8.705263157894738 15.789473684210526 1.5737982370014245 -0.5091934518542334 +9.492105263157894 15.789473684210526 2.5965352149740877 -0.8872935436291087 +10.278947368421052 15.789473684210526 4.396442976927523 -1.5846813552565104 +11.06578947368421 15.789473684210526 7.641847218147426 -2.900511000911099 +11.852631578947369 15.789473684210526 13.636895652018262 -5.440430742044206 +12.639473684210527 15.789473684210526 24.990727638679232 -10.456112617979404 +13.426315789473684 15.789473684210526 47.01714578365367 -20.593295087179104 +14.213157894736842 15.789473684210526 90.7849993566197 -41.54934565280752 +15.0 15.789473684210526 179.87254795797668 -85.86075456160873 +0.05 16.57894736842105 0.0008954684255083735 -0.009086135446008118 +0.5 16.57894736842105 0.009018685726437166 -0.009282424822437815 +0.8368421052631579 16.57894736842105 0.015291244964034034 -0.009646942267911577 +1.5 16.57894736842105 0.028657712501625774 -0.010962863889162568 +1.6236842105263158 16.57894736842105 0.031367075296160146 -0.011306006500792245 +2.4105263157894736 16.57894736842105 0.051012114292527395 -0.01436911537973599 +2.5 16.57894736842105 0.053576687347818695 -0.014829183059301405 +3.1973684210526314 16.57894736842105 0.07678577995455556 -0.019431696436432153 +3.5 16.57894736842105 0.08907195836805373 -0.022108663085421972 +3.9842105263157896 16.57894736842105 0.11246083951694469 -0.027530431405287104 +4.771052631578947 16.57894736842105 0.16389423398357664 -0.04043783501985364 +5.557894736842106 16.57894736842105 0.2404951429781756 -0.06121825341084383 +6.344736842105263 16.57894736842105 0.35782236639920484 -0.09515896875790893 +7.131578947368421 16.57894736842105 0.542132709444283 -0.1516309622316792 +7.918421052631579 16.57894736842105 0.8388059301055208 -0.24745671840404665 +8.705263157894738 16.57894736842105 1.3275694160106244 -0.41350245014085324 +9.492105263157894 16.57894736842105 2.1521667487827765 -0.7074796062539903 +10.278947368421052 16.57894736842105 3.5763583348606076 -1.2391674337442191 +11.06578947368421 16.57894736842105 6.094059325901349 -2.222090760410992 +11.852631578947369 16.57894736842105 10.651180576930221 -4.080033334354377 +12.639473684210527 16.57894736842105 19.09489958695614 -7.669899059981991 +13.426315789473684 16.57894736842105 35.1102168384886 -14.762613815704174 +14.213157894736842 16.57894736842105 66.21568010911355 -29.07995631007373 +15.0 16.57894736842105 128.01591647192475 -58.62407851231476 +0.05 17.36842105263158 0.0008364276093899223 -0.00848221236924663 +0.5 17.36842105263158 0.008421446681223951 -0.008657438424681644 +0.8368421052631579 17.36842105263158 0.014270571931430392 -0.008983180907135763 +1.5 17.36842105263158 0.02669316509273144 -0.010152826916279977 +1.6236842105263158 17.36842105263158 0.029202866657530507 -0.01045748627590907 +2.4105263157894736 17.36842105263158 0.04730464876157911 -0.013168353840859187 +2.5 17.36842105263158 0.04965605127100989 -0.013574202383238919 +3.1973684210526314 17.36842105263158 0.07081739331312602 -0.01761552374082955 +3.5 17.36842105263158 0.08194291326706035 -0.019956958525462427 +3.9842105263157896 17.36842105263158 0.10299905470827929 -0.02466665936844323 +4.771052631578947 17.36842105263158 0.1488440487323019 -0.03579300293322653 +5.557894736842106 17.36842105263158 0.21628424664394094 -0.053479807819885754 +6.344736842105263 17.36842105263158 0.31821916911383774 -0.08200827164837854 +7.131578947368421 17.36842105263158 0.4761324387585894 -0.12881092757596316 +7.918421052631579 17.36842105263158 0.7265820152290904 -0.20701365631182292 +8.705263157894738 17.36842105263158 1.1329129050472064 -0.3403656258806995 +9.492105263157894 17.36842105263158 1.807147048218457 -0.5724291887798639 +10.278947368421052 17.36842105263158 2.9514872083299752 -0.9846960463800187 +11.06578947368421 17.36842105263158 4.937921754949182 -1.7329850903214685 +11.852631578947369 17.36842105263158 8.46572291741952 -3.1201898457161628 +12.639473684210527 17.36842105263158 14.871261040555494 -5.745180758569281 +13.426315789473684 17.36842105263158 26.774787283496117 -10.823078761525958 +14.213157894736842 17.36842105263158 49.392966796920796 -20.853971754672106 +15.0 17.36842105263158 93.34811831265688 -41.09958825968762 +0.05 18.157894736842106 0.0007836039310756705 -0.007940739296834151 +0.5 18.157894736842106 0.007887357057085636 -0.008099542227043482 +0.8368421052631579 18.157894736842106 0.013358630055822098 -0.00839089100174107 +1.5 18.157894736842106 0.024943031197927956 -0.009437865944582723 +1.6236842105263158 18.157894736842106 0.02727587652213995 -0.009709896820067185 +2.4105263157894736 18.157894736842106 0.044023649433185964 -0.012122492932300884 +2.5 18.157894736842106 0.046188071221358776 -0.012482753167636435 +3.1973684210526314 18.157894736842106 0.06557538562050035 -0.016054663431785878 +3.5 18.157894736842106 0.07570402145087426 -0.018113460816429804 +3.9842105263157896 18.157894736842106 0.09476811376888711 -0.0222382816524567 +4.771052631578947 18.157894736842106 0.1358994382821802 -0.031898016232529254 +5.557894736842106 18.157894736842106 0.19569643076442755 -0.04709232631800145 +6.344736842105263 18.157894736842106 0.28497368061927963 -0.07130079814623971 +7.131578947368421 18.157894736842106 0.42152461774599137 -0.11048598503099795 +7.918421052631579 18.157894736842106 0.6351349854639462 -0.1750728167278211 +8.705263157894738 18.157894736842106 0.9767155728105807 -0.2835473567268476 +9.492105263157894 18.157894736842106 1.5348513008871747 -0.46937079616673405 +10.278947368421052 18.157894736842106 2.466969542758461 -0.7941414995936317 +11.06578947368421 18.157894736842106 4.057981033910274 -1.3732511145335944 +11.852631578947369 18.157894736842106 6.834019495009115 -2.4273520632801593 +12.639473684210527 18.157894736842106 11.781474252403564 -4.385220542955928 +13.426315789473684 18.157894736842106 20.80034544134105 -8.099212880909224 +14.213157894736842 18.157894736842106 37.595864127664036 -15.287711250480811 +15.0 18.157894736842106 69.55865864322226 -29.49333060769577 +0.05 18.94736842105263 0.0007361174136727883 -0.007456584100706878 +0.5 18.94736842105263 0.007407442720125345 -0.0075986763420962024 +0.8368421052631579 18.94736842105263 0.01253987372829643 -0.007860387057534666 +1.5 18.94736842105263 0.023376012087969737 -0.00880122356536319 +1.6236842105263158 18.94736842105263 0.025551980161015765 -0.009045327570448793 +2.4105263157894736 18.94736842105263 0.04110347920213849 -0.011203614898767388 +2.5 18.94736842105263 0.04310559589674073 -0.01152554130632316 +3.1973684210526314 18.94736842105263 0.06094400455620836 -0.014698521955953368 +3.5 18.94736842105263 0.07020971654313711 -0.01652078446278089 +3.9842105263157896 18.94736842105263 0.08755960926148955 -0.020155911205706853 +4.771052631578947 18.94736842105263 0.12466883121226466 -0.02860524297014387 +5.557894736842106 18.94736842105263 0.1780395437569714 -0.04175692835834711 +6.344736842105263 18.94736842105263 0.2568160382824192 -0.062468417744049494 +7.131578947368421 18.94736842105263 0.37584596951551774 -0.09559708942980688 +7.918421052631579 18.94736842105263 0.5596747723689148 -0.14950526573489475 +8.705263157894738 18.94736842105263 0.849724421868849 -0.23877637383397937 +9.492105263157894 18.94736842105263 1.3170332816771873 -0.38947540444820833 +10.278947368421052 18.94736842105263 2.085947755406089 -0.648803326689401 +11.06578947368421 18.94736842105263 3.3774820873816664 -1.103845316185561 +11.852631578947369 18.94736842105263 5.593695292723991 -1.9180973842247089 +12.639473684210527 18.94736842105263 9.476875688299685 -3.404146922787535 +13.426315789473684 18.94736842105263 16.427880424949723 -6.171719092164083 +14.213157894736842 18.94736842105263 29.13269910683986 -11.426614797079262 +15.0 18.94736842105263 52.85148195803469 -21.612460867362504 +0.05 19.736842105263158 0.0006932422105457152 -0.007017843993782561 +0.5 19.736842105263158 0.006974316593111866 -0.007146796665297263 +0.8368421052631579 19.736842105263158 0.01180146304553848 -0.007384473199849693 +1.5 19.736842105263158 0.021966355890965152 -0.008232494430474281 +1.6236842105263158 19.736842105263158 0.02400203969720924 -0.008452586826613212 +2.4105263157894736 19.736842105263158 0.03849221195572245 -0.01039314762278606 +2.5 19.736842105263158 0.040348948682943725 -0.010680310927616786 +3.1973684210526314 19.736842105263158 0.05682911674543752 -0.01351685176724009 +3.5 19.736842105263158 0.06534083428465215 -0.015138462909799087 +3.9842105263157896 19.736842105263158 0.08120518700444222 -0.01835989822936251 +4.771052631578947 19.736842105263158 0.11486489100293233 -0.025794241881785352 +5.557894736842106 19.736842105263158 0.16277807676832784 -0.03726038126349024 +6.344736842105263 19.736842105263158 0.2327502077697192 -0.05512131901053577 +7.131578947368421 19.736842105263158 0.33727345199405284 -0.08338589541182152 +7.918421052631579 19.736842105263158 0.4968092997466 -0.12878550674699796 +8.705263157894738 19.736842105263158 0.7454123493106594 -0.20304040520887257 +9.492105263157894 19.736842105263158 1.1405147910158444 -0.3266931715414788 +10.278947368421052 19.736842105263158 1.781758658820789 -0.5362991613453156 +11.06578947368421 19.736842105263158 2.8430754446157516 -0.8987470136856884 +11.852631578947369 19.736842105263158 4.636357951255684 -1.536910747543975 +12.639473684210527 19.736842105263158 7.728388174250893 -2.682855179795313 +13.426315789473684 19.736842105263158 13.169224756271909 -4.780527093727928 +14.213157894736842 19.736842105263158 22.943313480678675 -8.694011172190846 +15.0 19.736842105263158 40.86232591992561 -16.137504430895135 +0.05 20.526315789473685 0.00065437495918614 -0.006621869139532301 +0.5 20.526315789473685 0.006581824819549789 -0.006738436513696397 +0.8368421052631579 20.526315789473685 0.011132790263565258 -0.006953723989971284 +1.5 20.526315789473685 0.02069285777549912 -0.0077228075009234025 +1.6236842105263158 20.526315789473685 0.02260269995069924 -0.007922391838744665 +2.4105263157894736 20.526315789473685 0.036144423647159835 -0.009674357190110562 +2.5 20.526315789473685 0.03787333536244608 -0.009933089102049347 +3.1973684210526314 20.526315789473685 0.05315249963351295 -0.01247991211622793 +3.5 20.526315789473685 0.06100372588521825 -0.01392905519166874 +3.9842105263157896 20.526315789473685 0.07557118617454554 -0.016799006095407267 +4.771052631578947 20.526315789473685 0.10624572326201934 -0.0233796048388594 +5.557894736842106 20.526315789473685 0.1494914870150911 -0.033438882824814414 +6.344736842105263 20.526315789473685 0.21201141165635523 -0.04896390768832828 +7.131578947368421 20.526315789473685 0.3044350563831702 -0.07325402049074538 +7.918421052631579 20.526315789473685 0.44392421963594364 -0.11183801024083652 +8.705263157894738 20.526315789473685 0.6587357565277182 -0.17417990024998506 +9.492105263157894 20.526315789473685 0.9960031437198917 -0.27662975844318105 +10.278947368421052 20.526315789473685 1.5359873490617841 -0.44811971435034653 +11.06578947368421 20.526315789473685 2.4179471128275916 -0.74028820128123 +11.852631578947369 20.526315789473685 3.8862999384600116 -1.2474846326052096 +12.639473684210527 20.526315789473685 6.380454510791597 -2.1436130285964548 +13.426315789473684 20.526315789473685 10.700297410361328 -3.7579573521529412 +14.213157894736842 20.526315789473685 18.334011248505373 -6.720680172765909 +15.0 20.526315789473685 32.09265265301199 -12.258717030354685 +0.05 21.31578947368421 0.0006190095745122339 -0.006260899305580558 +0.5 21.31578947368421 0.0062248195119306874 -0.006367556762166132 +0.8368421052631579 21.31578947368421 0.010524937027160845 -0.006563693018925879 +1.5 21.31578947368421 0.019537660128566652 -0.007263890526476846 +1.6236842105263158 21.31578947368421 0.021334103797478337 -0.007444197178446523 +2.4105263157894736 21.31578947368421 0.03402546322257848 -0.009032578582839042 +2.5 21.31578947368421 0.03563974785181828 -0.009267098174863415 +3.1973684210526314 21.31578947368421 0.049852916923709974 -0.011562734985325694 +3.5 21.31578947368421 0.057120671947770166 -0.012865615194783126 +3.9842105263157896 21.31578947368421 0.07054909131596677 -0.015434844994858839 +4.771052631578947 21.31578947368421 0.09862565003079277 -0.02129075634932271 +5.557894736842106 21.31578947368421 0.13785337941202463 -0.030168456789778236 +6.344736842105263 21.31578947368421 0.19402431357310795 -0.04373752256654503 +7.131578947368421 21.31578947368421 0.27624126285937406 -0.06477016353473135 +7.918421052631579 21.31578947368421 0.39902730227417027 -0.09782287787224175 +8.705263157894738 21.31578947368421 0.5860311321580253 -0.1506364578147404 +9.492105263157894 21.31578947368421 0.876260389562841 -0.23636884502337546 +10.278947368421052 21.31578947368421 1.3353130515693674 -0.37807550194566375 +11.06578947368421 21.31578947368421 2.07530735510785 -0.6163155899011941 +11.852631578947369 21.31578947368421 3.2909638765418148 -1.0239486903514463 +12.639473684210527 21.31578947368421 5.326602333228009 -1.7341527528883534 +13.426315789473684 21.31578947368421 8.801282178936107 -2.9945225277740315 +14.213157894736842 21.31578947368421 14.846701432761906 -5.271380183518845 +15.0 21.31578947368421 25.568899208547766 -9.458140822540239 +0.05 22.105263157894736 0.0005867196872854669 -0.005932387602697752 +0.5 22.105263157894736 0.005898972576030648 -0.006029238885618811 +0.8368421052631579 22.105263157894736 0.00997050783999851 -0.006208762937749162 +1.5 22.105263157894736 0.018486228434163043 -0.0068478159415879675 +1.6236842105263158 22.105263157894736 0.020179793488839724 -0.007012364102590385 +2.4105263157894736 22.105263157894736 0.032105261331245374 -0.008456743536121832 +2.5 22.105263157894736 0.03361739391993851 -0.008670187248878379 +3.1973684210526314 22.105263157894736 0.04687933105341965 -0.010749440827758732 +3.5 22.105263157894736 0.053629598639180515 -0.011924430828023868 +3.9842105263157896 22.105263157894736 0.06605170010445034 -0.014233854088525 +4.771052631578947 22.105263157894736 0.09185247222599137 -0.01947086741957917 +5.557894736842106 22.105263157894736 0.1275894727677663 -0.02734444156263157 +6.344736842105263 22.105263157894736 0.17831588617445654 -0.039288914533892036 +7.131578947368421 22.105263157894736 0.2518655212158439 -0.0576191086794238 +7.918421052631579 22.105263157894736 0.36064544606395427 -0.08614338187024419 +8.705263157894738 22.105263157894736 0.5245837022672954 -0.13120761066629222 +9.492105263157894 22.105263157894736 0.7762220680871784 -0.20356756896252534 +10.278947368421052 22.105263157894736 1.1696859112803826 -0.32172307206640516 +11.06578947368421 22.105263157894736 1.7963842068538558 -0.5179973108945527 +11.852631578947369 22.105263157894736 2.8125994356331474 -0.8493793453482883 +12.639473684210527 22.105263157894736 4.492235769665198 -1.4191697015968825 +13.426315789473684 22.105263157894736 7.31957243727591 -2.4153311737751 +14.213157894736842 22.105263157894736 12.1667586593933 -4.188338063714941 +15.0 22.105263157894736 20.638431858773373 -7.401087132353606 +0.05 22.894736842105264 0.0005571433760154797 -0.00563081099298231 +0.5 22.894736842105264 0.005600600473573943 -0.005720000653323541 +0.8368421052631579 22.894736842105264 0.009463092600629004 -0.0058846479317923965 +1.5 22.894736842105264 0.017525552982632475 -0.006469924484486257 +1.6236842105263158 22.894736842105264 0.01912578881480675 -0.006620771959201002 +2.4105263157894736 22.894736842105264 0.030358572140142284 -0.007938995088460726 +2.5 22.894736842105264 0.03177804013596141 -0.008132691420597636 +3.1973684210526314 22.894736842105264 0.044188064413799195 -0.010021620974936883 +3.5 22.894736842105264 0.05047661679429641 -0.01108553445544386 +3.9842105263157896 22.894736842105264 0.06200723642226532 -0.013173986901517035 +4.771052631578947 22.894736842105264 0.0858022498677678 -0.017876238570023684 +5.557894736842106 22.894736842105264 0.11850012731630923 -0.024902337103910113 +6.344736842105263 22.894736842105264 0.16451393376252582 -0.03546694620512763 +7.131578947368421 22.894736842105264 0.2306522298273177 -0.05152863493023182 +7.918421052631579 22.894736842105264 0.3275459049315859 -0.07629999669665523 +8.705263157894738 22.894736842105264 0.4721421929217716 -0.1150696444109661 +9.492105263157894 22.894736842105264 0.6919126432449925 -0.1766514324554198 +10.278947368421052 22.894736842105264 1.031774493588329 -0.27606634316952894 +11.06578947368421 22.894736842105264 1.5667130497800534 -0.4391705442840245 +11.852631578947369 22.894736842105264 2.424243318681567 -0.7113395660840234 +12.639473684210527 22.894736842105264 3.823721810705477 -1.1730692484424838 +13.426315789473684 22.894736842105264 6.14848004040535 -1.9701412432559309 +14.213157894736842 22.894736842105264 10.081575485958469 -3.3689361260415422 +15.0 22.894736842105264 16.85452420946202 -5.866791393104007 +0.05 23.68421052631579 0.0005299718395687587 -0.005354404971752865 +0.5 23.68421052631579 0.005326566284476935 -0.005436816825863126 +0.8368421052631579 23.68421052631579 0.008997312467120202 -0.005587486555854928 +1.5 23.68421052631579 0.016645513580163328 -0.006125073355194881 +1.6236842105263158 23.68421052631579 0.01816061546393352 -0.006263922470873627 +2.4105263157894736 23.68421052631579 0.028764074085834933 -0.0074714951151326455 +2.5 23.68421052631579 0.03010053213639197 -0.007647606967749839 +3.1973684210526314 23.68421052631579 0.04174265018394783 -0.009370452099058827 +3.5 23.68421052631579 0.047618091763823596 -0.010339284044395974 +3.9842105263157896 23.68421052631579 0.058353235715590916 -0.012230189952872022 +4.771052631578947 23.68421052631579 0.08037223680857349 -0.016472677070991263 +5.557894736842106 23.68421052631579 0.11040052798534002 -0.022765809897179225 +6.344736842105263 23.68421052631579 0.15232668634682686 -0.0321545530469905 +7.131578947368421 23.68421052631579 0.2120780717256833 -0.04632239452568544 +7.918421052631579 23.68421052631579 0.2988558091646519 -0.06797585815388472 +8.705263157894738 23.68421052631579 0.42715227834957037 -0.10152431765546929 +9.492105263157894 23.68421052631579 0.6202094943096005 -0.15428208858412726 +10.278947368421052 23.68421052631579 0.9157768368975009 -0.2385898596831807 +11.06578947368421 23.68421052631579 1.3762012084741453 -0.3753534388527912 +11.852631578947369 23.68421052631579 2.1057565506371865 -0.6009835751396857 +12.639473684210527 23.68421052631579 3.282349453893203 -0.9790539744073617 +13.426315789473684 23.68421052631579 5.212428770523336 -1.6236149408742593 +14.213157894736842 23.68421052631579 8.436521759695822 -2.7394966068672337 +15.0 23.68421052631579 13.915098152146053 -4.705764184416182 +0.05 24.473684210526315 0.0005049397672626362 -0.0051004909679694945 +0.5 24.473684210526315 0.005074169462972422 -0.005176203691982657 +0.8368421052631579 24.473684210526315 0.008568514610619105 -0.0053151506618317424 +1.5 24.473684210526315 0.015836528039113017 -0.005810716967133103 +1.6236842105263158 24.473684210526315 0.017273690480720937 -0.005937734958710735 +2.4105263157894736 24.473684210526315 0.027304694948735643 -0.007046327720154506 +2.5 24.473684210526315 0.028564722477312603 -0.007208288343811512 +3.1973684210526314 24.473684210526315 0.039513489625015566 -0.008785045345553834 +3.5 24.473684210526315 0.04501682476378509 -0.009667002578651798 +3.9842105263157896 24.473684210526315 0.05503946379991069 -0.011390670549681856 +4.771052631578947 24.473684210526315 0.07548324881399905 -0.015232030344087006 +5.557894736842106 24.473684210526315 0.10315584260359684 -0.020891841140673648 +6.344736842105263 24.473684210526315 0.14150641669802944 -0.02928129836944618 +7.131578947368421 24.473684210526315 0.1957214885064367 -0.04184154226721097 +7.918421052631579 24.473684210526315 0.27381249164924826 -0.060869991291912624 +8.705263157894738 24.473684210526315 0.38829274932465546 -0.09008599436902255 +9.492105263157894 24.473684210526315 0.5589255103923116 -0.1356237039120202 +10.278947368421052 24.473684210526315 0.8176853548061294 -0.20759784773112547 +11.06578947368421 24.473684210526315 1.2165610530894337 -0.32324151266895146 +11.852631578947369 24.473684210526315 1.841967575072364 -0.511825103475029 +12.639473684210527 24.473684210526315 2.8394258421836605 -0.8243302517383491 +13.426315789473684 24.473684210526315 4.456614891809398 -1.3505029187125601 +14.213157894736842 24.473684210526315 7.125616146792478 -2.250646929789521 +15.0 24.473684210526315 11.603445225184661 -3.816214224352231 +0.05 25.263157894736842 0.00048181842366674314 -0.004864759099127755 +0.5 25.263157894736842 0.004841107129920398 -0.004935170358525633 +0.8368421052631579 25.263157894736842 0.008172754321629637 -0.0050639059813913975 +1.5 25.263157894736842 0.015091028547550466 -0.005521126223290092 +1.6236842105263158 25.263157894736842 0.016456747355386443 -0.0056387411544351175 +2.4105263157894736 25.263157894736842 0.025963769908822515 -0.006660348450643822 +2.5 25.263157894736842 0.027154850546773572 -0.006809787379307049 +3.1973684210526314 25.263157894736842 0.03747461094818256 -0.00825628189097387 +3.5 25.263157894736842 0.042642931018690655 -0.009063535034796738 +3.9842105263157896 25.263157894736842 0.05202427908151981 -0.010635448288306922 +4.771052631578947 25.263157894736842 0.07105712080425665 -0.014126096218329559 +5.557894736842106 25.263157894736842 0.09664713442462251 -0.019241435165578528 +6.344736842105263 25.263157894736842 0.13184998068659204 -0.02676437867076158 +7.131578947368421 25.263157894736842 0.18124647648473077 -0.037951573250233364 +7.918421052631579 25.263157894736842 0.2518332108032889 -0.0547582360115768 +8.705263157894738 25.263157894736842 0.35444357595631043 -0.08035072028733165 +9.492105263157894 25.263157894736842 0.506116146637667 -0.1198795443501057 +10.278947368421052 25.263157894736842 0.7339376179240047 -0.18179622266253512 +11.06578947368421 25.263157894736842 1.0818352484087161 -0.2802454706834918 +11.852631578947369 25.263157894736842 1.6216709120569737 -0.43915182328143737 +12.639473684210527 25.263157894736842 2.4737722497283015 -0.6996150065672352 +13.426315789473684 25.263157894736842 3.8400058998273967 -1.133257866424784 +14.213157894736842 25.263157894736842 6.0687364908805135 -1.8664272632288244 +15.0 25.263157894736842 9.763884255133398 -3.125732878315442 +0.05 26.05263157894737 0.0004604097778852859 -0.004647380734780913 +0.5 26.05263157894737 0.004625353792018123 -0.004712441846962487 +0.8368421052631579 26.05263157894737 0.00780652535367259 -0.004831826485211032 +1.5 26.05263157894737 0.014402246101032265 -0.005255162227571501 +1.6236842105263158 26.05263157894737 0.01570225406330191 -0.005363510487578426 +2.4105263157894736 26.05263157894737 0.0247293215082348 -0.0063068150548247544 +2.5 26.05263157894737 0.025857576819381962 -0.0064450128623550875 +3.1973684210526314 26.05263157894737 0.0356039876297597 -0.007775658372259726 +3.5 26.05263157894737 0.04046912534295574 -0.008516654432283837 +3.9842105263157896 26.05263157894737 0.049271737810025286 -0.00995753696200644 +4.771052631578947 26.05263157894737 0.06704119037349889 -0.013141617312898401 +5.557894736842106 26.05263157894737 0.09077253665937386 -0.01777439717753601 +6.344736842105263 26.05263157894737 0.12319858705148391 -0.024556923159718004 +7.131578947368421 26.05263157894737 0.16837228563052783 -0.034555194940497386 +7.918421052631579 26.05263157894737 0.23244769363443893 -0.049480489597330604 +8.705263157894738 26.05263157894737 0.3248435811734114 -0.07202456287528755 +9.492105263157894 26.05263157894737 0.46026358037559195 -0.10653689240316784 +10.278947368421052 26.05263157894737 0.662019006889101 -0.16012752345350667 +11.06578947368421 26.05263157894737 0.9672413482931925 -0.24452246334362185 +11.852631578947369 26.05263157894737 1.4363285201989178 -0.3793680698957248 +12.639473684210527 26.05263157894737 2.1692566439828016 -0.5981587760400692 +13.426315789473684 26.05263157894737 3.3325335180767275 -0.9585571690329748 +14.213157894736842 26.05263157894737 5.209066743857143 -1.5608943843960639 +15.0 26.05263157894737 8.285703211957467 -2.5835440757489585 +0.05 26.842105263157894 0.0004405411807885385 -0.004445940236249982 +0.5 26.842105263157894 0.004425167768715958 -0.004506133953489432 +0.8368421052631579 26.842105263157894 0.00746689454034104 -0.004616560762911825 +1.5 26.842105263157894 0.013764217021104056 -0.0050094406311340855 +1.6236842105263158 26.842105263157894 0.015003560289334346 -0.005110266020310056 +2.4105263157894736 26.842105263157894 0.02358915221831296 -0.00598347203552825 +2.5 26.842105263157894 0.024659398071832397 -0.006110809053908816 +3.1973684210526314 26.842105263157894 0.0338834646775021 -0.007338378083330451 +3.5 26.842105263157894 0.03847152883105839 -0.00802206682311351 +3.9842105263157896 26.842105263157894 0.04675148845929632 -0.009344158759814439 +4.771052631578947 26.842105263157894 0.06337910933609722 -0.012254976276217211 +5.557894736842106 26.842105263157894 0.08545320680191479 -0.016472882351365997 +6.344736842105263 26.842105263157894 0.1154149463509747 -0.02260237330914972 +7.131578947368421 26.842105263157894 0.1568684342417099 -0.031584851954405156 +7.918421052631579 26.842105263157894 0.21525083930083155 -0.044894587231728976 +8.705263157894738 26.842105263157894 0.29881570948142994 -0.0648505022490693 +9.492105263157894 26.842105263157894 0.42032827422615077 -0.0951428470138403 +10.278947368421052 26.842105263157894 0.5998498670886613 -0.1417637038721228 +11.06578947368421 26.842105263157894 0.8691348376837205 -0.2145514125639766 +11.852631578947369 26.842105263157894 1.279317247916555 -0.3297961688158551 +12.639473684210527 26.842105263157894 1.9139584963911733 -0.5149570539861872 +13.426315789473684 26.842105263157894 2.9110561181793906 -0.8166691450621733 +14.213157894736842 26.842105263157894 4.503547456626725 -1.3159612798207168 +15.0 26.842105263157894 7.085424917892163 -2.1544344387389684 +0.05 27.63157894736842 0.00042206189118505354 -0.004258281491265551 +0.5 27.63157894736842 0.004239024466971999 -0.004314678229739304 +0.8368421052631579 27.63157894736842 0.007151175219396086 -0.004417569770749997 +1.5 27.63157894736842 0.013172040825299113 -0.00478281745772999 +1.6236842105263158 27.63157894736842 0.014355265118214226 -0.0048767581082844935 +2.4105263157894736 27.63157894736842 0.02253389717169933 -0.005686124800726885 +2.5 27.63157894736842 0.023551314323261855 -0.00580450592764339 +3.1973684210526314 27.63157894736842 0.03229609860843198 -0.006940322097212255 +3.5 27.63157894736842 0.03663250669274162 -0.00757029854482953 +3.9842105263157896 27.63157894736842 0.04443662710919616 -0.008787916953885987 +4.771052631578947 27.63157894736842 0.06003580822554542 -0.01145949037504041 +5.557894736842106 27.63157894736842 0.08061939751496491 -0.015308207760969854 +6.344736842105263 27.63157894736842 0.1083853316763269 -0.020868459241377018 +7.131578947368421 27.63157894736842 0.146551445896942 -0.028968825045722884 +7.918421052631579 27.63157894736842 0.1999349316719778 -0.040888192790275224 +8.705263157894738 27.63157894736842 0.27581168275951934 -0.05862063744347348 +9.492105263157894 27.63157894736842 0.38533248935235864 -0.0853607604752115 +10.278947368421052 27.63157894736842 0.5458185736664958 -0.12616662692652364 +11.06578947368421 27.63157894736842 0.784543340980754 -0.1893409946516254 +11.852631578947369 27.63157894736842 1.1451535861933118 -0.2884515410424467 +12.639473684210527 27.63157894736842 1.6980306365973448 -0.4461822149240904 +13.426315789473684 27.63157894736842 2.558662945677238 -0.7007547682590666 +14.213157894736842 27.63157894736842 3.9191300217924154 -1.1177779299635107 +15.0 27.63157894736842 6.102842025815073 -1.8103289309611537 +0.05 28.42105263157895 0.00040483947655308057 -0.004083334230350588 +0.5 28.42105263157895 0.004065567220735519 -0.00413585941068865 +0.8368421052631579 28.42105263157895 0.006857097175488259 -0.004232098768259329 +1.5 28.42105263157895 0.012621072984170083 -0.004572438980336766 +1.6236842105263158 28.42105263157895 0.01375234600905341 -0.004659307825642305 +2.4105263157894736 28.42105263157895 0.02155476957363757 -0.005412861287090211 +2.5 28.42105263157895 0.022523504149534673 -0.005522456659777208 +3.1973684210526314 28.42105263157895 0.030828834044290494 -0.0065755673755331645 +3.5 28.42105263157895 0.03493480646895975 -0.007157701979087682 +3.9842105263157896 28.42105263157895 0.04230451214404033 -0.008282303523339545 +4.771052631578947 28.42105263157895 0.056971552596723325 -0.010740966836871124 +5.557894736842106 28.42105263157895 0.07621108407239731 -0.014264250145528828 +6.344736842105263 28.42105263157895 0.10201481238967566 -0.019326523610648495 +7.131578947368421 28.42105263157895 0.13726273192371116 -0.026657727985304653 +7.918421052631579 28.42105263157895 0.1862447405701858 -0.03737571308116737 +8.705263157894738 28.42105263157895 0.2553741072266109 -0.05321413378305539 +9.492105263157894 28.42105263157895 0.35446454107944303 -0.07689563261121675 +10.278947368421052 28.42105263157895 0.4986110965962023 -0.11279367583236324 +11.06578947368421 28.42105263157895 0.7113472656881416 -0.16788199972720327 +11.852631578947369 28.42105263157895 1.0299963232432283 -0.25364362135133006 +12.639473684210527 28.42105263157895 1.5143524186273116 -0.3888661115280537 +13.426315789473684 28.42105263157895 2.2616115166693147 -0.6051391676998468 +14.213157894736842 28.42105263157895 3.4316979721171923 -0.9558602108562707 +15.0 28.42105263157895 5.291626822977882 -1.5328358752036848 +0.05 29.210526315789473 0.000388757679517727 -0.0039201989243858575 +0.5 29.210526315789473 0.003903635180476895 -0.003969639781547456 +0.8368421052631579 29.210526315789473 0.006582617337988797 -0.004059262231821483 +1.5 29.210526315789473 0.01210741740378411 -0.004377076625688967 +1.6236842105263158 29.210526315789473 0.01319031221775589 -0.004458662138202314 +2.4105263157894736 29.210526315789473 0.0206447598921417 -0.005160118223028642 +2.5 29.210526315789473 0.021568297889779306 -0.0052625406276718125 +3.1973684210526314 29.210526315789473 0.029468688380913335 -0.006240621693259147 +3.5 29.210526315789473 0.03336355264021221 -0.0067805793796749765 +3.9842105263157896 29.210526315789473 0.04033577723695207 -0.0078219820959394 +4.771052631578947 29.210526315789473 0.05415324269375932 -0.010087662223209636 +5.557894736842106 29.210526315789473 0.07218130286789215 -0.013321458870712587 +6.344736842105263 29.210526315789473 0.0962181080928826 -0.017946511238821947 +7.131578947368421 29.210526315789473 0.12886539094065458 -0.024602833285307192 +7.918421052631579 29.210526315789473 0.17395149896259465 -0.0342694988059682 +8.705263157894738 29.210526315789473 0.23718233574272032 -0.04846956478716117 +9.492105263157894 29.210526315789473 0.3271537419995997 -0.069559642270348 +10.278947368421052 29.210526315789473 0.45715065557867535 -0.10127556795037108 +11.06578947368421 29.210526315789473 0.6475072857448944 -0.14960302060833683 +11.852631578947369 29.210526315789473 0.9304898244839024 -0.22415320102716366 +12.639473684210527 29.210526315789473 1.3570948576466961 -0.3407596138269055 +13.426315789473684 29.210526315789473 2.0095812748347814 -0.5256515699195755 +14.213157894736842 29.210526315789473 3.022303349408319 -0.8226683241970487 +15.0 29.210526315789473 4.616814725018867 -1.3068178845188165 +0.05 30.0 0.00037371321165078974 -0.003767644510680941 +0.5 30.0 0.003752172896788798 -0.0038137128642303742 +0.8368421052631579 30.0 0.006325991955397038 -0.0038978479301029463 +1.5 30.0 0.01162770781350448 -0.004195407797731633 +1.6236842105263158 30.0 0.0126656083448209 -0.0042711812958708215 +2.4105263157894736 30.0 0.019796550991729413 -0.004926349650071033 +2.5 30.0 0.020678143118648167 -0.005021675447950133 +3.1973684210526314 30.0 0.028205034276012352 -0.005932569137949413 +3.5 30.0 0.031905383534102646 -0.006435083509370783 +3.9842105263157896 30.0 0.03851419402575485 -0.0073992351809649155 +4.771052631578947 30.0 0.051557007234173266 -0.009495588857010474 +5.557894736842106 30.0 0.06848706207820293 -0.012472239187258335 +6.344736842105263 30.0 0.09093544969816426 -0.016706096994690115 +7.131578947368421 30.0 0.12124961207124589 -0.02277193399910646 +7.918421052631579 30.0 0.16286479642082535 -0.03152808378095451 +8.705263157894738 30.0 0.22086685086583988 -0.044303799059199546 +9.492105263157894 30.0 0.3028732868343674 -0.06315607436334905 +10.278947368421052 30.0 0.4205307890399163 -0.09129784224000985 +11.06578947368421 30.0 0.5916325988431066 -0.1338700173006138 +11.852631578947369 30.0 0.8441029392670165 -0.19905001564075675 +12.639473684210527 30.0 1.2216755122471643 -0.3001668053519415 +13.426315789473684 30.0 1.7944613739188713 -0.4591070192776733 +14.213157894736842 30.0 2.675692917269824 -0.7123037503297706 +15.0 30.0 4.0514307281349184 -1.121005048517917 diff --git a/test/enzyme_test.jl b/test/enzyme_test.jl index fcc2675..63f3abb 100644 --- a/test/enzyme_test.jl +++ b/test/enzyme_test.jl @@ -1,12 +1,16 @@ - -# This is just a placeholder at the moment to try and make CI trigger the error. -# This should hit the special method created in the EnzymeRules extension/weakdep. - using EnzymeCore, Enzyme import Bessels.BesselFunctions: besselkx_levin -dbesselk_dv(v, x) = autodiff(Forward, _v->besselkx_levin(_v, x, Val(20)), - Duplicated, Duplicated(v, 1.0))[2] +dbesselkx_dv(v, x) = autodiff(Forward, _v->besselkx_levin(_v, x, Val(30)), + Duplicated, Duplicated(v, 1.0))[2] -@test isapprox(dbesselk_dv(1.5, 13.1), 0.0410929, atol=1e-5) +dbesselkx_dx(v, x) = autodiff(Forward, _x->besselkx_levin(v, _x, Val(30)), + Duplicated, Duplicated(x, 1.0))[2] +for line in eachline("data/besselk/enzyme/besselkx_levin_enzyme_tests.csv") + (v, x, dv, dx) = parse.(Float64, split(line)) + test_dv = dbesselkx_dv(v, x) + test_dx = dbesselkx_dx(v, x) + @test isapprox(dv, test_dv, rtol=1e-4) + @test isapprox(dx, test_dx, rtol=1e-4) +end diff --git a/test/runtests.jl b/test/runtests.jl index db22996..2c1ec05 100644 --- a/test/runtests.jl +++ b/test/runtests.jl @@ -2,12 +2,12 @@ using Bessels using Test import SpecialFunctions -@time @testset "besseli" begin include("besseli_test.jl") end -@time @testset "besselk" begin include("besselk_test.jl") end -@time @testset "besselj" begin include("besselj_test.jl") end -@time @testset "bessely" begin include("bessely_test.jl") end -@time @testset "hankel" begin include("hankel_test.jl") end -@time @testset "gamma" begin include("gamma_test.jl") end -@time @testset "airy" begin include("airy_test.jl") end -@time @testset "sphericalbessel" begin include("sphericalbessel_test.jl") end +#@time @testset "besseli" begin include("besseli_test.jl") end +#@time @testset "besselk" begin include("besselk_test.jl") end +#@time @testset "besselj" begin include("besselj_test.jl") end +#@time @testset "bessely" begin include("bessely_test.jl") end +#@time @testset "hankel" begin include("hankel_test.jl") end +#@time @testset "gamma" begin include("gamma_test.jl") end +#@time @testset "airy" begin include("airy_test.jl") end +#@time @testset "sphericalbessel" begin include("sphericalbessel_test.jl") end @time @testset "enzyme autodiff" begin include("enzyme_test.jl") end From 71334c4343e49c5a1901e5ad1be4d9dda640ddb4 Mon Sep 17 00:00:00 2001 From: Chris Geoga Date: Tue, 9 May 2023 11:46:00 -0400 Subject: [PATCH 08/17] Oops, for got to un-comment the other tests. --- test/runtests.jl | 16 ++++++++-------- 1 file changed, 8 insertions(+), 8 deletions(-) diff --git a/test/runtests.jl b/test/runtests.jl index 2c1ec05..db22996 100644 --- a/test/runtests.jl +++ b/test/runtests.jl @@ -2,12 +2,12 @@ using Bessels using Test import SpecialFunctions -#@time @testset "besseli" begin include("besseli_test.jl") end -#@time @testset "besselk" begin include("besselk_test.jl") end -#@time @testset "besselj" begin include("besselj_test.jl") end -#@time @testset "bessely" begin include("bessely_test.jl") end -#@time @testset "hankel" begin include("hankel_test.jl") end -#@time @testset "gamma" begin include("gamma_test.jl") end -#@time @testset "airy" begin include("airy_test.jl") end -#@time @testset "sphericalbessel" begin include("sphericalbessel_test.jl") end +@time @testset "besseli" begin include("besseli_test.jl") end +@time @testset "besselk" begin include("besselk_test.jl") end +@time @testset "besselj" begin include("besselj_test.jl") end +@time @testset "bessely" begin include("bessely_test.jl") end +@time @testset "hankel" begin include("hankel_test.jl") end +@time @testset "gamma" begin include("gamma_test.jl") end +@time @testset "airy" begin include("airy_test.jl") end +@time @testset "sphericalbessel" begin include("sphericalbessel_test.jl") end @time @testset "enzyme autodiff" begin include("enzyme_test.jl") end From afcb93489c7dad331539af578a57c219e3ab2d1b Mon Sep 17 00:00:00 2001 From: Chris Geoga Date: Tue, 9 May 2023 12:17:21 -0400 Subject: [PATCH 09/17] Switch to simple fd for enzyme test. Rtol=1e-15 works for all but six (v,x) pairs. Setting to 5e-14 for now, which seems pretty acceptable considering the possibility of edge cases in the Arb function too. --- .../enzyme/besselkx_levin_enzyme_tests.csv | 960 +++++++++--------- test/enzyme_test.jl | 4 +- 2 files changed, 482 insertions(+), 482 deletions(-) diff --git a/test/data/besselk/enzyme/besselkx_levin_enzyme_tests.csv b/test/data/besselk/enzyme/besselkx_levin_enzyme_tests.csv index 4cdc411..2812537 100644 --- a/test/data/besselk/enzyme/besselkx_levin_enzyme_tests.csv +++ b/test/data/besselk/enzyme/besselkx_levin_enzyme_tests.csv @@ -1,480 +1,480 @@ -0.05 15.0 0.0010367551275860735 -0.010535444408076956 -0.5 15.0 0.01044909892351376 -0.010786835957018058 -0.8368421052631579 15.0 0.017739559706085498 -0.011253628340044715 -1.5 15.0 0.03339462144739582 -0.01294415009942762 -1.6236842105263158 15.0 0.03659268674207609 -0.013387189370521184 -2.4105263157894736 15.0 0.060052964058098846 -0.017373888634610528 -2.5 15.0 0.06315118081428582 -0.017978120147579314 -3.1973684210526314 15.0 0.09153794927366792 -0.024071625383024185 -3.5 15.0 0.10679821808137248 -0.027662120748450387 -3.9842105263157896 15.0 0.1362234972398303 -0.03500862169499547 -4.771052631578947 15.0 0.20238301287210217 -0.05287136846744076 -5.557894736842106 15.0 0.30374695994764694 -0.08239500819064587 -6.344736842105263 15.0 0.4636380330289938 -0.13208397871236013 -7.131578947368421 15.0 0.7228219860633713 -0.21742481203781536 -7.918421052631579 15.0 1.1539993886345965 -0.36729242248277877 -8.705263157894738 15.0 1.8899224810295472 -0.6365147679147345 -9.492105263157894 15.0 3.178648748030142 -1.1317368209637213 -10.278947368421052 15.0 5.4925111585656206 -2.064343411173481 -11.06578947368421 15.0 9.755179997903623 -3.862762417806891 -11.852631578947369 15.0 17.809093701857034 -7.414088150681344 -12.639473684210527 15.0 33.42266223980047 -14.596901632387276 -13.426315789473684 15.0 64.44838838752855 -29.471567616652393 -14.213157894736842 15.0 127.68673796168575 -61.01851035096151 -15.0 15.0 259.78299436259425 -129.48144922416114 -0.05 15.789473684210526 0.0009618017615888898 -0.009765556568247924 -0.5 15.789473684210526 0.009690053231145072 -0.009988315614370639 -0.8368421052631579 15.789473684210526 0.016439733325060442 -0.010399134668978188 -1.5 15.789473684210526 0.03087566747671769 -0.011886188577293443 -1.6236842105263158 15.789473684210526 0.03381256980329193 -0.0122738720493344 -2.4105263157894736 15.789473684210526 0.05522747048187236 -0.015757523702255714 -2.5 15.789473684210526 0.058038666965933475 -0.016282009216853487 -3.1973684210526314 15.789473684210526 0.08362918041182572 -0.02155695614180951 -3.5 15.789473684210526 0.09727720780815732 -0.024645549572446225 -3.9842105263157896 15.789473684210526 0.12342043461990765 -0.03092882160405939 -4.771052631578947 15.789473684210526 0.18152361114544602 -0.046042595633940836 -5.557894736842106 15.789473684210526 0.2692415361359061 -0.07068159896315275 -6.344736842105263 15.789473684210526 0.4055175223058745 -0.11150082248932047 -7.131578947368421 15.789473684210526 0.622835091385322 -0.18046448445800095 -7.918421052631579 15.789473684210526 0.9782051119542616 -0.29944048626010167 -8.705263157894738 15.789473684210526 1.5737982370014245 -0.5091934518542334 -9.492105263157894 15.789473684210526 2.5965352149740877 -0.8872935436291087 -10.278947368421052 15.789473684210526 4.396442976927523 -1.5846813552565104 -11.06578947368421 15.789473684210526 7.641847218147426 -2.900511000911099 -11.852631578947369 15.789473684210526 13.636895652018262 -5.440430742044206 -12.639473684210527 15.789473684210526 24.990727638679232 -10.456112617979404 -13.426315789473684 15.789473684210526 47.01714578365367 -20.593295087179104 -14.213157894736842 15.789473684210526 90.7849993566197 -41.54934565280752 -15.0 15.789473684210526 179.87254795797668 -85.86075456160873 -0.05 16.57894736842105 0.0008954684255083735 -0.009086135446008118 -0.5 16.57894736842105 0.009018685726437166 -0.009282424822437815 -0.8368421052631579 16.57894736842105 0.015291244964034034 -0.009646942267911577 -1.5 16.57894736842105 0.028657712501625774 -0.010962863889162568 -1.6236842105263158 16.57894736842105 0.031367075296160146 -0.011306006500792245 -2.4105263157894736 16.57894736842105 0.051012114292527395 -0.01436911537973599 -2.5 16.57894736842105 0.053576687347818695 -0.014829183059301405 -3.1973684210526314 16.57894736842105 0.07678577995455556 -0.019431696436432153 -3.5 16.57894736842105 0.08907195836805373 -0.022108663085421972 -3.9842105263157896 16.57894736842105 0.11246083951694469 -0.027530431405287104 -4.771052631578947 16.57894736842105 0.16389423398357664 -0.04043783501985364 -5.557894736842106 16.57894736842105 0.2404951429781756 -0.06121825341084383 -6.344736842105263 16.57894736842105 0.35782236639920484 -0.09515896875790893 -7.131578947368421 16.57894736842105 0.542132709444283 -0.1516309622316792 -7.918421052631579 16.57894736842105 0.8388059301055208 -0.24745671840404665 -8.705263157894738 16.57894736842105 1.3275694160106244 -0.41350245014085324 -9.492105263157894 16.57894736842105 2.1521667487827765 -0.7074796062539903 -10.278947368421052 16.57894736842105 3.5763583348606076 -1.2391674337442191 -11.06578947368421 16.57894736842105 6.094059325901349 -2.222090760410992 -11.852631578947369 16.57894736842105 10.651180576930221 -4.080033334354377 -12.639473684210527 16.57894736842105 19.09489958695614 -7.669899059981991 -13.426315789473684 16.57894736842105 35.1102168384886 -14.762613815704174 -14.213157894736842 16.57894736842105 66.21568010911355 -29.07995631007373 -15.0 16.57894736842105 128.01591647192475 -58.62407851231476 -0.05 17.36842105263158 0.0008364276093899223 -0.00848221236924663 -0.5 17.36842105263158 0.008421446681223951 -0.008657438424681644 -0.8368421052631579 17.36842105263158 0.014270571931430392 -0.008983180907135763 -1.5 17.36842105263158 0.02669316509273144 -0.010152826916279977 -1.6236842105263158 17.36842105263158 0.029202866657530507 -0.01045748627590907 -2.4105263157894736 17.36842105263158 0.04730464876157911 -0.013168353840859187 -2.5 17.36842105263158 0.04965605127100989 -0.013574202383238919 -3.1973684210526314 17.36842105263158 0.07081739331312602 -0.01761552374082955 -3.5 17.36842105263158 0.08194291326706035 -0.019956958525462427 -3.9842105263157896 17.36842105263158 0.10299905470827929 -0.02466665936844323 -4.771052631578947 17.36842105263158 0.1488440487323019 -0.03579300293322653 -5.557894736842106 17.36842105263158 0.21628424664394094 -0.053479807819885754 -6.344736842105263 17.36842105263158 0.31821916911383774 -0.08200827164837854 -7.131578947368421 17.36842105263158 0.4761324387585894 -0.12881092757596316 -7.918421052631579 17.36842105263158 0.7265820152290904 -0.20701365631182292 -8.705263157894738 17.36842105263158 1.1329129050472064 -0.3403656258806995 -9.492105263157894 17.36842105263158 1.807147048218457 -0.5724291887798639 -10.278947368421052 17.36842105263158 2.9514872083299752 -0.9846960463800187 -11.06578947368421 17.36842105263158 4.937921754949182 -1.7329850903214685 -11.852631578947369 17.36842105263158 8.46572291741952 -3.1201898457161628 -12.639473684210527 17.36842105263158 14.871261040555494 -5.745180758569281 -13.426315789473684 17.36842105263158 26.774787283496117 -10.823078761525958 -14.213157894736842 17.36842105263158 49.392966796920796 -20.853971754672106 -15.0 17.36842105263158 93.34811831265688 -41.09958825968762 -0.05 18.157894736842106 0.0007836039310756705 -0.007940739296834151 -0.5 18.157894736842106 0.007887357057085636 -0.008099542227043482 -0.8368421052631579 18.157894736842106 0.013358630055822098 -0.00839089100174107 -1.5 18.157894736842106 0.024943031197927956 -0.009437865944582723 -1.6236842105263158 18.157894736842106 0.02727587652213995 -0.009709896820067185 -2.4105263157894736 18.157894736842106 0.044023649433185964 -0.012122492932300884 -2.5 18.157894736842106 0.046188071221358776 -0.012482753167636435 -3.1973684210526314 18.157894736842106 0.06557538562050035 -0.016054663431785878 -3.5 18.157894736842106 0.07570402145087426 -0.018113460816429804 -3.9842105263157896 18.157894736842106 0.09476811376888711 -0.0222382816524567 -4.771052631578947 18.157894736842106 0.1358994382821802 -0.031898016232529254 -5.557894736842106 18.157894736842106 0.19569643076442755 -0.04709232631800145 -6.344736842105263 18.157894736842106 0.28497368061927963 -0.07130079814623971 -7.131578947368421 18.157894736842106 0.42152461774599137 -0.11048598503099795 -7.918421052631579 18.157894736842106 0.6351349854639462 -0.1750728167278211 -8.705263157894738 18.157894736842106 0.9767155728105807 -0.2835473567268476 -9.492105263157894 18.157894736842106 1.5348513008871747 -0.46937079616673405 -10.278947368421052 18.157894736842106 2.466969542758461 -0.7941414995936317 -11.06578947368421 18.157894736842106 4.057981033910274 -1.3732511145335944 -11.852631578947369 18.157894736842106 6.834019495009115 -2.4273520632801593 -12.639473684210527 18.157894736842106 11.781474252403564 -4.385220542955928 -13.426315789473684 18.157894736842106 20.80034544134105 -8.099212880909224 -14.213157894736842 18.157894736842106 37.595864127664036 -15.287711250480811 -15.0 18.157894736842106 69.55865864322226 -29.49333060769577 -0.05 18.94736842105263 0.0007361174136727883 -0.007456584100706878 -0.5 18.94736842105263 0.007407442720125345 -0.0075986763420962024 -0.8368421052631579 18.94736842105263 0.01253987372829643 -0.007860387057534666 -1.5 18.94736842105263 0.023376012087969737 -0.00880122356536319 -1.6236842105263158 18.94736842105263 0.025551980161015765 -0.009045327570448793 -2.4105263157894736 18.94736842105263 0.04110347920213849 -0.011203614898767388 -2.5 18.94736842105263 0.04310559589674073 -0.01152554130632316 -3.1973684210526314 18.94736842105263 0.06094400455620836 -0.014698521955953368 -3.5 18.94736842105263 0.07020971654313711 -0.01652078446278089 -3.9842105263157896 18.94736842105263 0.08755960926148955 -0.020155911205706853 -4.771052631578947 18.94736842105263 0.12466883121226466 -0.02860524297014387 -5.557894736842106 18.94736842105263 0.1780395437569714 -0.04175692835834711 -6.344736842105263 18.94736842105263 0.2568160382824192 -0.062468417744049494 -7.131578947368421 18.94736842105263 0.37584596951551774 -0.09559708942980688 -7.918421052631579 18.94736842105263 0.5596747723689148 -0.14950526573489475 -8.705263157894738 18.94736842105263 0.849724421868849 -0.23877637383397937 -9.492105263157894 18.94736842105263 1.3170332816771873 -0.38947540444820833 -10.278947368421052 18.94736842105263 2.085947755406089 -0.648803326689401 -11.06578947368421 18.94736842105263 3.3774820873816664 -1.103845316185561 -11.852631578947369 18.94736842105263 5.593695292723991 -1.9180973842247089 -12.639473684210527 18.94736842105263 9.476875688299685 -3.404146922787535 -13.426315789473684 18.94736842105263 16.427880424949723 -6.171719092164083 -14.213157894736842 18.94736842105263 29.13269910683986 -11.426614797079262 -15.0 18.94736842105263 52.85148195803469 -21.612460867362504 -0.05 19.736842105263158 0.0006932422105457152 -0.007017843993782561 -0.5 19.736842105263158 0.006974316593111866 -0.007146796665297263 -0.8368421052631579 19.736842105263158 0.01180146304553848 -0.007384473199849693 -1.5 19.736842105263158 0.021966355890965152 -0.008232494430474281 -1.6236842105263158 19.736842105263158 0.02400203969720924 -0.008452586826613212 -2.4105263157894736 19.736842105263158 0.03849221195572245 -0.01039314762278606 -2.5 19.736842105263158 0.040348948682943725 -0.010680310927616786 -3.1973684210526314 19.736842105263158 0.05682911674543752 -0.01351685176724009 -3.5 19.736842105263158 0.06534083428465215 -0.015138462909799087 -3.9842105263157896 19.736842105263158 0.08120518700444222 -0.01835989822936251 -4.771052631578947 19.736842105263158 0.11486489100293233 -0.025794241881785352 -5.557894736842106 19.736842105263158 0.16277807676832784 -0.03726038126349024 -6.344736842105263 19.736842105263158 0.2327502077697192 -0.05512131901053577 -7.131578947368421 19.736842105263158 0.33727345199405284 -0.08338589541182152 -7.918421052631579 19.736842105263158 0.4968092997466 -0.12878550674699796 -8.705263157894738 19.736842105263158 0.7454123493106594 -0.20304040520887257 -9.492105263157894 19.736842105263158 1.1405147910158444 -0.3266931715414788 -10.278947368421052 19.736842105263158 1.781758658820789 -0.5362991613453156 -11.06578947368421 19.736842105263158 2.8430754446157516 -0.8987470136856884 -11.852631578947369 19.736842105263158 4.636357951255684 -1.536910747543975 -12.639473684210527 19.736842105263158 7.728388174250893 -2.682855179795313 -13.426315789473684 19.736842105263158 13.169224756271909 -4.780527093727928 -14.213157894736842 19.736842105263158 22.943313480678675 -8.694011172190846 -15.0 19.736842105263158 40.86232591992561 -16.137504430895135 -0.05 20.526315789473685 0.00065437495918614 -0.006621869139532301 -0.5 20.526315789473685 0.006581824819549789 -0.006738436513696397 -0.8368421052631579 20.526315789473685 0.011132790263565258 -0.006953723989971284 -1.5 20.526315789473685 0.02069285777549912 -0.0077228075009234025 -1.6236842105263158 20.526315789473685 0.02260269995069924 -0.007922391838744665 -2.4105263157894736 20.526315789473685 0.036144423647159835 -0.009674357190110562 -2.5 20.526315789473685 0.03787333536244608 -0.009933089102049347 -3.1973684210526314 20.526315789473685 0.05315249963351295 -0.01247991211622793 -3.5 20.526315789473685 0.06100372588521825 -0.01392905519166874 -3.9842105263157896 20.526315789473685 0.07557118617454554 -0.016799006095407267 -4.771052631578947 20.526315789473685 0.10624572326201934 -0.0233796048388594 -5.557894736842106 20.526315789473685 0.1494914870150911 -0.033438882824814414 -6.344736842105263 20.526315789473685 0.21201141165635523 -0.04896390768832828 -7.131578947368421 20.526315789473685 0.3044350563831702 -0.07325402049074538 -7.918421052631579 20.526315789473685 0.44392421963594364 -0.11183801024083652 -8.705263157894738 20.526315789473685 0.6587357565277182 -0.17417990024998506 -9.492105263157894 20.526315789473685 0.9960031437198917 -0.27662975844318105 -10.278947368421052 20.526315789473685 1.5359873490617841 -0.44811971435034653 -11.06578947368421 20.526315789473685 2.4179471128275916 -0.74028820128123 -11.852631578947369 20.526315789473685 3.8862999384600116 -1.2474846326052096 -12.639473684210527 20.526315789473685 6.380454510791597 -2.1436130285964548 -13.426315789473684 20.526315789473685 10.700297410361328 -3.7579573521529412 -14.213157894736842 20.526315789473685 18.334011248505373 -6.720680172765909 -15.0 20.526315789473685 32.09265265301199 -12.258717030354685 -0.05 21.31578947368421 0.0006190095745122339 -0.006260899305580558 -0.5 21.31578947368421 0.0062248195119306874 -0.006367556762166132 -0.8368421052631579 21.31578947368421 0.010524937027160845 -0.006563693018925879 -1.5 21.31578947368421 0.019537660128566652 -0.007263890526476846 -1.6236842105263158 21.31578947368421 0.021334103797478337 -0.007444197178446523 -2.4105263157894736 21.31578947368421 0.03402546322257848 -0.009032578582839042 -2.5 21.31578947368421 0.03563974785181828 -0.009267098174863415 -3.1973684210526314 21.31578947368421 0.049852916923709974 -0.011562734985325694 -3.5 21.31578947368421 0.057120671947770166 -0.012865615194783126 -3.9842105263157896 21.31578947368421 0.07054909131596677 -0.015434844994858839 -4.771052631578947 21.31578947368421 0.09862565003079277 -0.02129075634932271 -5.557894736842106 21.31578947368421 0.13785337941202463 -0.030168456789778236 -6.344736842105263 21.31578947368421 0.19402431357310795 -0.04373752256654503 -7.131578947368421 21.31578947368421 0.27624126285937406 -0.06477016353473135 -7.918421052631579 21.31578947368421 0.39902730227417027 -0.09782287787224175 -8.705263157894738 21.31578947368421 0.5860311321580253 -0.1506364578147404 -9.492105263157894 21.31578947368421 0.876260389562841 -0.23636884502337546 -10.278947368421052 21.31578947368421 1.3353130515693674 -0.37807550194566375 -11.06578947368421 21.31578947368421 2.07530735510785 -0.6163155899011941 -11.852631578947369 21.31578947368421 3.2909638765418148 -1.0239486903514463 -12.639473684210527 21.31578947368421 5.326602333228009 -1.7341527528883534 -13.426315789473684 21.31578947368421 8.801282178936107 -2.9945225277740315 -14.213157894736842 21.31578947368421 14.846701432761906 -5.271380183518845 -15.0 21.31578947368421 25.568899208547766 -9.458140822540239 -0.05 22.105263157894736 0.0005867196872854669 -0.005932387602697752 -0.5 22.105263157894736 0.005898972576030648 -0.006029238885618811 -0.8368421052631579 22.105263157894736 0.00997050783999851 -0.006208762937749162 -1.5 22.105263157894736 0.018486228434163043 -0.0068478159415879675 -1.6236842105263158 22.105263157894736 0.020179793488839724 -0.007012364102590385 -2.4105263157894736 22.105263157894736 0.032105261331245374 -0.008456743536121832 -2.5 22.105263157894736 0.03361739391993851 -0.008670187248878379 -3.1973684210526314 22.105263157894736 0.04687933105341965 -0.010749440827758732 -3.5 22.105263157894736 0.053629598639180515 -0.011924430828023868 -3.9842105263157896 22.105263157894736 0.06605170010445034 -0.014233854088525 -4.771052631578947 22.105263157894736 0.09185247222599137 -0.01947086741957917 -5.557894736842106 22.105263157894736 0.1275894727677663 -0.02734444156263157 -6.344736842105263 22.105263157894736 0.17831588617445654 -0.039288914533892036 -7.131578947368421 22.105263157894736 0.2518655212158439 -0.0576191086794238 -7.918421052631579 22.105263157894736 0.36064544606395427 -0.08614338187024419 -8.705263157894738 22.105263157894736 0.5245837022672954 -0.13120761066629222 -9.492105263157894 22.105263157894736 0.7762220680871784 -0.20356756896252534 -10.278947368421052 22.105263157894736 1.1696859112803826 -0.32172307206640516 -11.06578947368421 22.105263157894736 1.7963842068538558 -0.5179973108945527 -11.852631578947369 22.105263157894736 2.8125994356331474 -0.8493793453482883 -12.639473684210527 22.105263157894736 4.492235769665198 -1.4191697015968825 -13.426315789473684 22.105263157894736 7.31957243727591 -2.4153311737751 -14.213157894736842 22.105263157894736 12.1667586593933 -4.188338063714941 -15.0 22.105263157894736 20.638431858773373 -7.401087132353606 -0.05 22.894736842105264 0.0005571433760154797 -0.00563081099298231 -0.5 22.894736842105264 0.005600600473573943 -0.005720000653323541 -0.8368421052631579 22.894736842105264 0.009463092600629004 -0.0058846479317923965 -1.5 22.894736842105264 0.017525552982632475 -0.006469924484486257 -1.6236842105263158 22.894736842105264 0.01912578881480675 -0.006620771959201002 -2.4105263157894736 22.894736842105264 0.030358572140142284 -0.007938995088460726 -2.5 22.894736842105264 0.03177804013596141 -0.008132691420597636 -3.1973684210526314 22.894736842105264 0.044188064413799195 -0.010021620974936883 -3.5 22.894736842105264 0.05047661679429641 -0.01108553445544386 -3.9842105263157896 22.894736842105264 0.06200723642226532 -0.013173986901517035 -4.771052631578947 22.894736842105264 0.0858022498677678 -0.017876238570023684 -5.557894736842106 22.894736842105264 0.11850012731630923 -0.024902337103910113 -6.344736842105263 22.894736842105264 0.16451393376252582 -0.03546694620512763 -7.131578947368421 22.894736842105264 0.2306522298273177 -0.05152863493023182 -7.918421052631579 22.894736842105264 0.3275459049315859 -0.07629999669665523 -8.705263157894738 22.894736842105264 0.4721421929217716 -0.1150696444109661 -9.492105263157894 22.894736842105264 0.6919126432449925 -0.1766514324554198 -10.278947368421052 22.894736842105264 1.031774493588329 -0.27606634316952894 -11.06578947368421 22.894736842105264 1.5667130497800534 -0.4391705442840245 -11.852631578947369 22.894736842105264 2.424243318681567 -0.7113395660840234 -12.639473684210527 22.894736842105264 3.823721810705477 -1.1730692484424838 -13.426315789473684 22.894736842105264 6.14848004040535 -1.9701412432559309 -14.213157894736842 22.894736842105264 10.081575485958469 -3.3689361260415422 -15.0 22.894736842105264 16.85452420946202 -5.866791393104007 -0.05 23.68421052631579 0.0005299718395687587 -0.005354404971752865 -0.5 23.68421052631579 0.005326566284476935 -0.005436816825863126 -0.8368421052631579 23.68421052631579 0.008997312467120202 -0.005587486555854928 -1.5 23.68421052631579 0.016645513580163328 -0.006125073355194881 -1.6236842105263158 23.68421052631579 0.01816061546393352 -0.006263922470873627 -2.4105263157894736 23.68421052631579 0.028764074085834933 -0.0074714951151326455 -2.5 23.68421052631579 0.03010053213639197 -0.007647606967749839 -3.1973684210526314 23.68421052631579 0.04174265018394783 -0.009370452099058827 -3.5 23.68421052631579 0.047618091763823596 -0.010339284044395974 -3.9842105263157896 23.68421052631579 0.058353235715590916 -0.012230189952872022 -4.771052631578947 23.68421052631579 0.08037223680857349 -0.016472677070991263 -5.557894736842106 23.68421052631579 0.11040052798534002 -0.022765809897179225 -6.344736842105263 23.68421052631579 0.15232668634682686 -0.0321545530469905 -7.131578947368421 23.68421052631579 0.2120780717256833 -0.04632239452568544 -7.918421052631579 23.68421052631579 0.2988558091646519 -0.06797585815388472 -8.705263157894738 23.68421052631579 0.42715227834957037 -0.10152431765546929 -9.492105263157894 23.68421052631579 0.6202094943096005 -0.15428208858412726 -10.278947368421052 23.68421052631579 0.9157768368975009 -0.2385898596831807 -11.06578947368421 23.68421052631579 1.3762012084741453 -0.3753534388527912 -11.852631578947369 23.68421052631579 2.1057565506371865 -0.6009835751396857 -12.639473684210527 23.68421052631579 3.282349453893203 -0.9790539744073617 -13.426315789473684 23.68421052631579 5.212428770523336 -1.6236149408742593 -14.213157894736842 23.68421052631579 8.436521759695822 -2.7394966068672337 -15.0 23.68421052631579 13.915098152146053 -4.705764184416182 -0.05 24.473684210526315 0.0005049397672626362 -0.0051004909679694945 -0.5 24.473684210526315 0.005074169462972422 -0.005176203691982657 -0.8368421052631579 24.473684210526315 0.008568514610619105 -0.0053151506618317424 -1.5 24.473684210526315 0.015836528039113017 -0.005810716967133103 -1.6236842105263158 24.473684210526315 0.017273690480720937 -0.005937734958710735 -2.4105263157894736 24.473684210526315 0.027304694948735643 -0.007046327720154506 -2.5 24.473684210526315 0.028564722477312603 -0.007208288343811512 -3.1973684210526314 24.473684210526315 0.039513489625015566 -0.008785045345553834 -3.5 24.473684210526315 0.04501682476378509 -0.009667002578651798 -3.9842105263157896 24.473684210526315 0.05503946379991069 -0.011390670549681856 -4.771052631578947 24.473684210526315 0.07548324881399905 -0.015232030344087006 -5.557894736842106 24.473684210526315 0.10315584260359684 -0.020891841140673648 -6.344736842105263 24.473684210526315 0.14150641669802944 -0.02928129836944618 -7.131578947368421 24.473684210526315 0.1957214885064367 -0.04184154226721097 -7.918421052631579 24.473684210526315 0.27381249164924826 -0.060869991291912624 -8.705263157894738 24.473684210526315 0.38829274932465546 -0.09008599436902255 -9.492105263157894 24.473684210526315 0.5589255103923116 -0.1356237039120202 -10.278947368421052 24.473684210526315 0.8176853548061294 -0.20759784773112547 -11.06578947368421 24.473684210526315 1.2165610530894337 -0.32324151266895146 -11.852631578947369 24.473684210526315 1.841967575072364 -0.511825103475029 -12.639473684210527 24.473684210526315 2.8394258421836605 -0.8243302517383491 -13.426315789473684 24.473684210526315 4.456614891809398 -1.3505029187125601 -14.213157894736842 24.473684210526315 7.125616146792478 -2.250646929789521 -15.0 24.473684210526315 11.603445225184661 -3.816214224352231 -0.05 25.263157894736842 0.00048181842366674314 -0.004864759099127755 -0.5 25.263157894736842 0.004841107129920398 -0.004935170358525633 -0.8368421052631579 25.263157894736842 0.008172754321629637 -0.0050639059813913975 -1.5 25.263157894736842 0.015091028547550466 -0.005521126223290092 -1.6236842105263158 25.263157894736842 0.016456747355386443 -0.0056387411544351175 -2.4105263157894736 25.263157894736842 0.025963769908822515 -0.006660348450643822 -2.5 25.263157894736842 0.027154850546773572 -0.006809787379307049 -3.1973684210526314 25.263157894736842 0.03747461094818256 -0.00825628189097387 -3.5 25.263157894736842 0.042642931018690655 -0.009063535034796738 -3.9842105263157896 25.263157894736842 0.05202427908151981 -0.010635448288306922 -4.771052631578947 25.263157894736842 0.07105712080425665 -0.014126096218329559 -5.557894736842106 25.263157894736842 0.09664713442462251 -0.019241435165578528 -6.344736842105263 25.263157894736842 0.13184998068659204 -0.02676437867076158 -7.131578947368421 25.263157894736842 0.18124647648473077 -0.037951573250233364 -7.918421052631579 25.263157894736842 0.2518332108032889 -0.0547582360115768 -8.705263157894738 25.263157894736842 0.35444357595631043 -0.08035072028733165 -9.492105263157894 25.263157894736842 0.506116146637667 -0.1198795443501057 -10.278947368421052 25.263157894736842 0.7339376179240047 -0.18179622266253512 -11.06578947368421 25.263157894736842 1.0818352484087161 -0.2802454706834918 -11.852631578947369 25.263157894736842 1.6216709120569737 -0.43915182328143737 -12.639473684210527 25.263157894736842 2.4737722497283015 -0.6996150065672352 -13.426315789473684 25.263157894736842 3.8400058998273967 -1.133257866424784 -14.213157894736842 25.263157894736842 6.0687364908805135 -1.8664272632288244 -15.0 25.263157894736842 9.763884255133398 -3.125732878315442 -0.05 26.05263157894737 0.0004604097778852859 -0.004647380734780913 -0.5 26.05263157894737 0.004625353792018123 -0.004712441846962487 -0.8368421052631579 26.05263157894737 0.00780652535367259 -0.004831826485211032 -1.5 26.05263157894737 0.014402246101032265 -0.005255162227571501 -1.6236842105263158 26.05263157894737 0.01570225406330191 -0.005363510487578426 -2.4105263157894736 26.05263157894737 0.0247293215082348 -0.0063068150548247544 -2.5 26.05263157894737 0.025857576819381962 -0.0064450128623550875 -3.1973684210526314 26.05263157894737 0.0356039876297597 -0.007775658372259726 -3.5 26.05263157894737 0.04046912534295574 -0.008516654432283837 -3.9842105263157896 26.05263157894737 0.049271737810025286 -0.00995753696200644 -4.771052631578947 26.05263157894737 0.06704119037349889 -0.013141617312898401 -5.557894736842106 26.05263157894737 0.09077253665937386 -0.01777439717753601 -6.344736842105263 26.05263157894737 0.12319858705148391 -0.024556923159718004 -7.131578947368421 26.05263157894737 0.16837228563052783 -0.034555194940497386 -7.918421052631579 26.05263157894737 0.23244769363443893 -0.049480489597330604 -8.705263157894738 26.05263157894737 0.3248435811734114 -0.07202456287528755 -9.492105263157894 26.05263157894737 0.46026358037559195 -0.10653689240316784 -10.278947368421052 26.05263157894737 0.662019006889101 -0.16012752345350667 -11.06578947368421 26.05263157894737 0.9672413482931925 -0.24452246334362185 -11.852631578947369 26.05263157894737 1.4363285201989178 -0.3793680698957248 -12.639473684210527 26.05263157894737 2.1692566439828016 -0.5981587760400692 -13.426315789473684 26.05263157894737 3.3325335180767275 -0.9585571690329748 -14.213157894736842 26.05263157894737 5.209066743857143 -1.5608943843960639 -15.0 26.05263157894737 8.285703211957467 -2.5835440757489585 -0.05 26.842105263157894 0.0004405411807885385 -0.004445940236249982 -0.5 26.842105263157894 0.004425167768715958 -0.004506133953489432 -0.8368421052631579 26.842105263157894 0.00746689454034104 -0.004616560762911825 -1.5 26.842105263157894 0.013764217021104056 -0.0050094406311340855 -1.6236842105263158 26.842105263157894 0.015003560289334346 -0.005110266020310056 -2.4105263157894736 26.842105263157894 0.02358915221831296 -0.00598347203552825 -2.5 26.842105263157894 0.024659398071832397 -0.006110809053908816 -3.1973684210526314 26.842105263157894 0.0338834646775021 -0.007338378083330451 -3.5 26.842105263157894 0.03847152883105839 -0.00802206682311351 -3.9842105263157896 26.842105263157894 0.04675148845929632 -0.009344158759814439 -4.771052631578947 26.842105263157894 0.06337910933609722 -0.012254976276217211 -5.557894736842106 26.842105263157894 0.08545320680191479 -0.016472882351365997 -6.344736842105263 26.842105263157894 0.1154149463509747 -0.02260237330914972 -7.131578947368421 26.842105263157894 0.1568684342417099 -0.031584851954405156 -7.918421052631579 26.842105263157894 0.21525083930083155 -0.044894587231728976 -8.705263157894738 26.842105263157894 0.29881570948142994 -0.0648505022490693 -9.492105263157894 26.842105263157894 0.42032827422615077 -0.0951428470138403 -10.278947368421052 26.842105263157894 0.5998498670886613 -0.1417637038721228 -11.06578947368421 26.842105263157894 0.8691348376837205 -0.2145514125639766 -11.852631578947369 26.842105263157894 1.279317247916555 -0.3297961688158551 -12.639473684210527 26.842105263157894 1.9139584963911733 -0.5149570539861872 -13.426315789473684 26.842105263157894 2.9110561181793906 -0.8166691450621733 -14.213157894736842 26.842105263157894 4.503547456626725 -1.3159612798207168 -15.0 26.842105263157894 7.085424917892163 -2.1544344387389684 -0.05 27.63157894736842 0.00042206189118505354 -0.004258281491265551 -0.5 27.63157894736842 0.004239024466971999 -0.004314678229739304 -0.8368421052631579 27.63157894736842 0.007151175219396086 -0.004417569770749997 -1.5 27.63157894736842 0.013172040825299113 -0.00478281745772999 -1.6236842105263158 27.63157894736842 0.014355265118214226 -0.0048767581082844935 -2.4105263157894736 27.63157894736842 0.02253389717169933 -0.005686124800726885 -2.5 27.63157894736842 0.023551314323261855 -0.00580450592764339 -3.1973684210526314 27.63157894736842 0.03229609860843198 -0.006940322097212255 -3.5 27.63157894736842 0.03663250669274162 -0.00757029854482953 -3.9842105263157896 27.63157894736842 0.04443662710919616 -0.008787916953885987 -4.771052631578947 27.63157894736842 0.06003580822554542 -0.01145949037504041 -5.557894736842106 27.63157894736842 0.08061939751496491 -0.015308207760969854 -6.344736842105263 27.63157894736842 0.1083853316763269 -0.020868459241377018 -7.131578947368421 27.63157894736842 0.146551445896942 -0.028968825045722884 -7.918421052631579 27.63157894736842 0.1999349316719778 -0.040888192790275224 -8.705263157894738 27.63157894736842 0.27581168275951934 -0.05862063744347348 -9.492105263157894 27.63157894736842 0.38533248935235864 -0.0853607604752115 -10.278947368421052 27.63157894736842 0.5458185736664958 -0.12616662692652364 -11.06578947368421 27.63157894736842 0.784543340980754 -0.1893409946516254 -11.852631578947369 27.63157894736842 1.1451535861933118 -0.2884515410424467 -12.639473684210527 27.63157894736842 1.6980306365973448 -0.4461822149240904 -13.426315789473684 27.63157894736842 2.558662945677238 -0.7007547682590666 -14.213157894736842 27.63157894736842 3.9191300217924154 -1.1177779299635107 -15.0 27.63157894736842 6.102842025815073 -1.8103289309611537 -0.05 28.42105263157895 0.00040483947655308057 -0.004083334230350588 -0.5 28.42105263157895 0.004065567220735519 -0.00413585941068865 -0.8368421052631579 28.42105263157895 0.006857097175488259 -0.004232098768259329 -1.5 28.42105263157895 0.012621072984170083 -0.004572438980336766 -1.6236842105263158 28.42105263157895 0.01375234600905341 -0.004659307825642305 -2.4105263157894736 28.42105263157895 0.02155476957363757 -0.005412861287090211 -2.5 28.42105263157895 0.022523504149534673 -0.005522456659777208 -3.1973684210526314 28.42105263157895 0.030828834044290494 -0.0065755673755331645 -3.5 28.42105263157895 0.03493480646895975 -0.007157701979087682 -3.9842105263157896 28.42105263157895 0.04230451214404033 -0.008282303523339545 -4.771052631578947 28.42105263157895 0.056971552596723325 -0.010740966836871124 -5.557894736842106 28.42105263157895 0.07621108407239731 -0.014264250145528828 -6.344736842105263 28.42105263157895 0.10201481238967566 -0.019326523610648495 -7.131578947368421 28.42105263157895 0.13726273192371116 -0.026657727985304653 -7.918421052631579 28.42105263157895 0.1862447405701858 -0.03737571308116737 -8.705263157894738 28.42105263157895 0.2553741072266109 -0.05321413378305539 -9.492105263157894 28.42105263157895 0.35446454107944303 -0.07689563261121675 -10.278947368421052 28.42105263157895 0.4986110965962023 -0.11279367583236324 -11.06578947368421 28.42105263157895 0.7113472656881416 -0.16788199972720327 -11.852631578947369 28.42105263157895 1.0299963232432283 -0.25364362135133006 -12.639473684210527 28.42105263157895 1.5143524186273116 -0.3888661115280537 -13.426315789473684 28.42105263157895 2.2616115166693147 -0.6051391676998468 -14.213157894736842 28.42105263157895 3.4316979721171923 -0.9558602108562707 -15.0 28.42105263157895 5.291626822977882 -1.5328358752036848 -0.05 29.210526315789473 0.000388757679517727 -0.0039201989243858575 -0.5 29.210526315789473 0.003903635180476895 -0.003969639781547456 -0.8368421052631579 29.210526315789473 0.006582617337988797 -0.004059262231821483 -1.5 29.210526315789473 0.01210741740378411 -0.004377076625688967 -1.6236842105263158 29.210526315789473 0.01319031221775589 -0.004458662138202314 -2.4105263157894736 29.210526315789473 0.0206447598921417 -0.005160118223028642 -2.5 29.210526315789473 0.021568297889779306 -0.0052625406276718125 -3.1973684210526314 29.210526315789473 0.029468688380913335 -0.006240621693259147 -3.5 29.210526315789473 0.03336355264021221 -0.0067805793796749765 -3.9842105263157896 29.210526315789473 0.04033577723695207 -0.0078219820959394 -4.771052631578947 29.210526315789473 0.05415324269375932 -0.010087662223209636 -5.557894736842106 29.210526315789473 0.07218130286789215 -0.013321458870712587 -6.344736842105263 29.210526315789473 0.0962181080928826 -0.017946511238821947 -7.131578947368421 29.210526315789473 0.12886539094065458 -0.024602833285307192 -7.918421052631579 29.210526315789473 0.17395149896259465 -0.0342694988059682 -8.705263157894738 29.210526315789473 0.23718233574272032 -0.04846956478716117 -9.492105263157894 29.210526315789473 0.3271537419995997 -0.069559642270348 -10.278947368421052 29.210526315789473 0.45715065557867535 -0.10127556795037108 -11.06578947368421 29.210526315789473 0.6475072857448944 -0.14960302060833683 -11.852631578947369 29.210526315789473 0.9304898244839024 -0.22415320102716366 -12.639473684210527 29.210526315789473 1.3570948576466961 -0.3407596138269055 -13.426315789473684 29.210526315789473 2.0095812748347814 -0.5256515699195755 -14.213157894736842 29.210526315789473 3.022303349408319 -0.8226683241970487 -15.0 29.210526315789473 4.616814725018867 -1.3068178845188165 -0.05 30.0 0.00037371321165078974 -0.003767644510680941 -0.5 30.0 0.003752172896788798 -0.0038137128642303742 -0.8368421052631579 30.0 0.006325991955397038 -0.0038978479301029463 -1.5 30.0 0.01162770781350448 -0.004195407797731633 -1.6236842105263158 30.0 0.0126656083448209 -0.0042711812958708215 -2.4105263157894736 30.0 0.019796550991729413 -0.004926349650071033 -2.5 30.0 0.020678143118648167 -0.005021675447950133 -3.1973684210526314 30.0 0.028205034276012352 -0.005932569137949413 -3.5 30.0 0.031905383534102646 -0.006435083509370783 -3.9842105263157896 30.0 0.03851419402575485 -0.0073992351809649155 -4.771052631578947 30.0 0.051557007234173266 -0.009495588857010474 -5.557894736842106 30.0 0.06848706207820293 -0.012472239187258335 -6.344736842105263 30.0 0.09093544969816426 -0.016706096994690115 -7.131578947368421 30.0 0.12124961207124589 -0.02277193399910646 -7.918421052631579 30.0 0.16286479642082535 -0.03152808378095451 -8.705263157894738 30.0 0.22086685086583988 -0.044303799059199546 -9.492105263157894 30.0 0.3028732868343674 -0.06315607436334905 -10.278947368421052 30.0 0.4205307890399163 -0.09129784224000985 -11.06578947368421 30.0 0.5916325988431066 -0.1338700173006138 -11.852631578947369 30.0 0.8441029392670165 -0.19905001564075675 -12.639473684210527 30.0 1.2216755122471643 -0.3001668053519415 -13.426315789473684 30.0 1.7944613739188713 -0.4591070192776733 -14.213157894736842 30.0 2.675692917269824 -0.7123037503297706 -15.0 30.0 4.0514307281349184 -1.121005048517917 +0.05 15.0 0.0010367550243248873 -0.010535067751984675 +0.5 15.0 0.010449098914014986 -0.01078681062530944 +0.8368421052631579 15.0 0.017739524427007214 -0.011253247948315986 +1.5 15.0 0.03339475018149044 -0.012944172750371328 +1.6236842105263158 15.0 0.03659256330986051 -0.013387134424641038 +2.4105263157894736 15.0 0.060053461161278285 -0.01737441813255621 +2.5 15.0 0.06315177428496668 -0.01797801770884907 +3.1973684210526314 15.0 0.09153796353729869 -0.02407158316293611 +3.5 15.0 0.10679732917798115 -0.027662176581349097 +3.9842105263157896 15.0 0.13622339317144114 -0.03500871760185614 +4.771052631578947 15.0 0.20238655469920025 -0.052870923772866806 +5.557894736842106 15.0 0.303740141100424 -0.0823982000399553 +6.344736842105263 15.0 0.46363797403353085 -0.13208326320796876 +7.131578947368421 15.0 0.7228249251400523 -0.21742629587213247 +7.918421052631579 15.0 1.1540273601164805 -0.36728954520731494 +8.705263157894738 15.0 1.8899943448433185 -0.6365448043443811 +9.492105263157894 15.0 3.178504728294413 -1.1317443315551583 +10.278947368421052 15.0 5.492508162447618 -2.064269307995206 +11.06578947368421 15.0 9.755248128546953 -3.8625840980194326 +11.852631578947369 15.0 17.810002387947012 -7.414119464537326 +12.639473684210527 15.0 33.42072274175519 -14.597025494943916 +13.426315789473684 15.0 64.44813807960297 -29.472816282994675 +14.213157894736842 15.0 127.6808854062801 -61.015296240751546 +15.0 15.0 259.7840931825677 -129.48031158793447 +0.05 15.789473684210526 0.0009618018465905567 -0.009765937006862346 +0.5 15.789473684210526 0.009690049078591587 -0.009987998356122674 +0.8368421052631579 15.789473684210526 0.016439702321741138 -0.010399087245865534 +1.5 15.789473684210526 0.03087564745420991 -0.011885718043785982 +1.6236842105263158 15.789473684210526 0.033812578142526924 -0.012274269652290662 +2.4105263157894736 15.789473684210526 0.05522752834932757 -0.015756999670943225 +2.5 15.789473684210526 0.05803870843626657 -0.016282101986872646 +3.1973684210526314 15.789473684210526 0.08362996107785413 -0.021556865877259922 +3.5 15.789473684210526 0.09727819962171295 -0.024645458079276383 +3.9842105263157896 15.789473684210526 0.12341929317453491 -0.0309278881378392 +4.771052631578947 15.789473684210526 0.18152309437490025 -0.046042278706613655 +5.557894736842106 15.789473684210526 0.269242314076902 -0.07067884069081469 +6.344736842105263 15.789473684210526 0.4055175837170212 -0.11150082793616066 +7.131578947368421 15.789473684210526 0.6228462752991946 -0.18046339407618 +7.918421052631579 15.789473684210526 0.9782251482109025 -0.29942851214966415 +8.705263157894738 15.789473684210526 1.5737927323383258 -0.509175273763553 +9.492105263157894 15.789473684210526 2.5965429158562605 -0.8873256034755608 +10.278947368421052 15.789473684210526 4.396248916630848 -1.584691961533538 +11.06578947368421 15.789473684210526 7.641414852807336 -2.9004050394268637 +11.852631578947369 15.789473684210526 13.637642516684103 -5.440214017236562 +12.639473684210527 15.789473684210526 24.990621970407172 -10.456534869297457 +13.426315789473684 15.789473684210526 47.01457037999508 -20.593151349288803 +14.213157894736842 15.789473684210526 90.78444985403569 -41.547853307692264 +15.0 15.789473684210526 179.88315762214583 -85.85659078228966 +0.05 16.57894736842105 0.0008954685031868151 -0.00908607229325002 +0.5 16.57894736842105 0.009018684024561801 -0.009283131737476141 +0.8368421052631579 16.57894736842105 0.015291246785994835 -0.009647650957644676 +1.5 16.57894736842105 0.02865782704110998 -0.010962936528067062 +1.6236842105263158 16.57894736842105 0.031367217953135726 -0.011305986788562851 +2.4105263157894736 16.57894736842105 0.05101207685287158 -0.014369025060696164 +2.5 16.57894736842105 0.053576652731264635 -0.014829153903236641 +3.1973684210526314 16.57894736842105 0.0767856692739055 -0.019430266588248434 +3.5 16.57894736842105 0.08907289751387044 -0.022108900519532956 +3.9842105263157896 16.57894736842105 0.11245980601364239 -0.02752822186202752 +4.771052631578947 16.57894736842105 0.16389076226830493 -0.04044041642998249 +5.557894736842106 16.57894736842105 0.2404953140391126 -0.06121755200775037 +6.344736842105263 16.57894736842105 0.3578233075498357 -0.09516025223304479 +7.131578947368421 16.57894736842105 0.5421442165915208 -0.15162982178870016 +7.918421052631579 16.57894736842105 0.8387898343603494 -0.2474608488181228 +8.705263157894738 16.57894736842105 1.327622487832204 -0.4135064618590878 +9.492105263157894 16.57894736842105 2.152262583047653 -0.707420224401091 +10.278947368421052 16.57894736842105 3.5763703979090002 -1.23907769746871 +11.06578947368421 16.57894736842105 6.094105219318691 -2.2220772927863934 +11.852631578947369 16.57894736842105 10.651093774932527 -4.080025876528624 +12.639473684210527 16.57894736842105 19.095020952706125 -7.67000320855779 +13.426315789473684 16.57894736842105 35.11227394292208 -14.761223925251796 +14.213157894736842 16.57894736842105 66.21271998516191 -29.079605014990833 +15.0 16.57894736842105 128.01703057724623 -58.63007772548719 +0.05 17.36842105263158 0.0008364276988424413 -0.008481609569931332 +0.5 17.36842105263158 0.00842144499790296 -0.00865743888546313 +0.8368421052631579 17.36842105263158 0.014270575368674524 -0.00898245451505397 +1.5 17.36842105263158 0.0266931816195864 -0.010152814692952218 +1.6236842105263158 17.36842105263158 0.029202738086453434 -0.010457482417872829 +2.4105263157894736 17.36842105263158 0.047304604101216424 -0.013168218663786169 +2.5 17.36842105263158 0.049655675411305854 -0.013574053282813612 +3.1973684210526314 17.36842105263158 0.07081673577560098 -0.017615353579991132 +3.5 17.36842105263158 0.08194364562001111 -0.01995562790075154 +3.9842105263157896 17.36842105263158 0.10299890042334789 -0.02466703659585353 +4.771052631578947 17.36842105263158 0.14884628124075377 -0.035792894486788934 +5.557894736842106 17.36842105263158 0.21628386753179543 -0.05348386484346484 +6.344736842105263 17.36842105263158 0.3182181990126527 -0.08200924216465101 +7.131578947368421 17.36842105263158 0.47614133505973727 -0.12879929276132335 +7.918421052631579 17.36842105263158 0.7265858808037061 -0.20701010118129073 +8.705263157894738 17.36842105263158 1.132902396340547 -0.3403620371401612 +9.492105263157894 17.36842105263158 1.8071431329523409 -0.5724266201859811 +10.278947368421052 17.36842105263158 2.9514726261410553 -0.9847677080349462 +11.06578947368421 17.36842105263158 4.937977040995968 -1.7330135483326807 +11.852631578947369 17.36842105263158 8.465334909262017 -3.1198752876502183 +12.639473684210527 17.36842105263158 14.872107686423767 -5.74562885598651 +13.426315789473684 17.36842105263158 26.77493741249963 -10.823865233258264 +14.213157894736842 17.36842105263158 49.39295945219402 -20.856048089185425 +15.0 17.36842105263158 93.34776368489538 -41.09877441560721 +0.05 18.157894736842106 0.0007836039205990373 -0.007941336135992466 +0.5 18.157894736842106 0.0078873554957884 -0.008099005645026731 +0.8368421052631579 18.157894736842106 0.013358628627561266 -0.008390262507555988 +1.5 18.157894736842106 0.024943043184086584 -0.009437102229857235 +1.6236842105263158 18.157894736842106 0.027275995971163307 -0.009709136731588568 +2.4105263157894736 18.157894736842106 0.04402365419505333 -0.01212168520099287 +2.5 18.157894736842106 0.046188531381735815 -0.01248175677794983 +3.1973684210526314 18.157894736842106 0.06557533291488388 -0.016053613839385913 +3.5 18.157894736842106 0.07570480231137645 -0.01811193657756306 +3.9842105263157896 18.157894736842106 0.09476878987906492 -0.022236861518794792 +4.771052631578947 18.157894736842106 0.13589937195597732 -0.031897879714991174 +5.557894736842106 18.157894736842106 0.19569664071895532 -0.04709146781702443 +6.344736842105263 18.157894736842106 0.2849794118570284 -0.07129497370014763 +7.131578947368421 18.157894736842106 0.42152271538909997 -0.11047846524522838 +7.918421052631579 18.157894736842106 0.6351192732453554 -0.17506091857813752 +8.705263157894738 18.157894736842106 0.9766767854712329 -0.2835438797293878 +9.492105263157894 18.157894736842106 1.534859375774044 -0.4693742515935451 +10.278947368421052 18.157894736842106 2.467069288740017 -0.7941271073244253 +11.06578947368421 18.157894736842106 4.0581534418994325 -1.373265721923606 +11.852631578947369 18.157894736842106 6.833637294669379 -2.4273532403085283 +12.639473684210527 18.157894736842106 11.782061832963285 -4.385641511876809 +13.426315789473684 18.157894736842106 20.799431597832626 -8.099278892051794 +14.213157894736842 18.157894736842106 37.59374266037649 -15.287842656026426 +15.0 18.157894736842106 69.55920824047091 -29.49096478706292 +0.05 18.94736842105263 0.0007361174070886092 -0.007456096079949244 +0.5 18.94736842105263 0.007407444466740315 -0.007598127783295703 +0.8368421052631579 18.94736842105263 0.012539853104345964 -0.007860339366565914 +1.5 18.94736842105263 0.023376015124409347 -0.008801164682317522 +1.6236842105263158 18.94736842105263 0.025552004102557414 -0.009045254342315 +2.4105263157894736 18.94736842105263 0.04110378018279115 -0.011203510878198965 +2.5 18.94736842105263 0.043105207193317036 -0.011524706550936365 +3.1973684210526314 18.94736842105263 0.06094404549079722 -0.014699650152215991 +3.5 18.94736842105263 0.07020963354170938 -0.01652097634526514 +3.9842105263157896 18.94736842105263 0.08755953769443692 -0.02015557393859916 +4.771052631578947 18.94736842105263 0.12467108228540681 -0.028603513958625048 +5.557894736842106 18.94736842105263 0.17803977883699393 -0.04175430977209193 +6.344736842105263 18.94736842105263 0.25681588521595944 -0.06246916470721465 +7.131578947368421 18.94736842105263 0.3758451272159583 -0.09559828551637488 +7.918421052631579 18.94736842105263 0.5596857490730777 -0.14949285673145263 +8.705263157894738 18.94736842105263 0.8497278797231034 -0.238773246473716 +9.492105263157894 18.94736842105263 1.3170276350207428 -0.38947891707887555 +10.278947368421052 18.94736842105263 2.085839243587362 -0.6488068003323374 +11.06578947368421 18.94736842105263 3.3775042636385932 -1.1038329166769656 +11.852631578947369 18.94736842105263 5.593723278449159 -1.9181079309001061 +12.639473684210527 18.94736842105263 9.477294379226311 -3.404398763458295 +13.426315789473684 18.94736842105263 16.427800858277635 -6.171728909288849 +14.213157894736842 18.94736842105263 29.132472232306586 -11.427657239988472 +15.0 18.94736842105263 52.84899005346882 -21.61026170774502 +0.05 19.736842105263158 0.0006932422843103742 -0.007018347144651731 +0.5 19.736842105263158 0.006974315863502267 -0.007146829850702864 +0.8368421052631579 19.736842105263158 0.011801443969107168 -0.007383895926482416 +1.5 19.736842105263158 0.021966368876393305 -0.0082331479880097 +1.6236842105263158 19.736842105263158 0.0240021231959649 -0.008453149946727433 +2.4105263157894736 19.736842105263158 0.038491868878694456 -0.010393051689247508 +2.5 19.736842105263158 0.04034894751860129 -0.010680984857407769 +3.1973684210526314 19.736842105263158 0.05682859793244545 -0.01351781166014191 +3.5 19.736842105263158 0.06534085467563862 -0.015138346126082735 +3.9842105263157896 19.736842105263158 0.08120461527741767 -0.018359595996760687 +4.771052631578947 19.736842105263158 0.11486476186048343 -0.02579392507172195 +5.557894736842106 19.736842105263158 0.16277770816045875 -0.03725736392006179 +6.344736842105263 19.736842105263158 0.2327457208057145 -0.0551260015958644 +7.131578947368421 19.736842105263158 0.3372789710932994 -0.08338010427797599 +7.918421052631579 19.736842105263158 0.49681898197580293 -0.12878701950574484 +8.705263157894738 19.736842105263158 0.7453766808088381 -0.2030377291590343 +9.492105263157894 19.736842105263158 1.1405575786440934 -0.32666577040687234 +10.278947368421052 19.736842105263158 1.7816827341310364 -0.5363494540926786 +11.06578947368421 19.736842105263158 2.84308708744328 -0.8987369652601022 +11.852631578947369 19.736842105263158 4.63632908520415 -1.5370455755057972 +12.639473684210527 19.736842105263158 7.728366457931344 -2.683073405144214 +13.426315789473684 19.736842105263158 13.169837816849125 -4.7805707216016105 +14.213157894736842 19.736842105263158 22.943510610061956 -8.694067194483097 +15.0 19.736842105263158 40.86020225634043 -16.137694067920915 +0.05 20.526315789473685 0.0006543748889245913 -0.006621826186075766 +0.5 20.526315789473685 0.0065818233036432595 -0.006738503119279088 +0.8368421052631579 20.526315789473685 0.01113277372086845 -0.006953675098516337 +1.5 20.526315789473685 0.020692842103650588 -0.007723361267481417 +1.6236842105263158 20.526315789473685 0.022602709851817754 -0.007922474607203437 +2.4105263157894736 20.526315789473685 0.0361444207775913 -0.009673692750547653 +2.5 20.526315789473685 0.037873318439665495 -0.009932978907678948 +3.1973684210526314 20.526315789473685 0.05315254242908035 -0.012479790268554065 +3.5 20.526315789473685 0.06100376017622138 -0.01392897123700731 +3.9842105263157896 20.526315789473685 0.0755705064172324 -0.016799118037548164 +4.771052631578947 20.526315789473685 0.10624548739485264 -0.023379562853636977 +5.557894736842106 20.526315789473685 0.14949184202580623 -0.03343676407514754 +6.344736842105263 20.526315789473685 0.2120120602670754 -0.04896068221668 +7.131578947368421 20.526315789473685 0.3044342724563555 -0.07324782776767133 +7.918421052631579 20.526315789473685 0.44392505906313956 -0.11183704741226988 +8.705263157894738 20.526315789473685 0.6587031210529125 -0.17417832214129936 +9.492105263157894 20.526315789473685 0.9959622445454742 -0.27665324620729986 +10.278947368421052 20.526315789473685 1.5360013234786178 -0.4481257288593995 +11.06578947368421 20.526315789473685 2.417828897505321 -0.74030091603591 +11.852631578947369 20.526315789473685 3.8863098175500417 -1.247360651166101 +12.639473684210527 20.526315789473685 6.380389536371569 -2.1437708085731617 +13.426315789473684 20.526315789473685 10.700884402934639 -3.7582343429658427 +14.213157894736842 20.526315789473685 18.33483444758378 -6.7206543435788735 +15.0 20.526315789473685 32.092830843567945 -12.258856548829591 +0.05 21.31578947368421 0.0006190095212039973 -0.006261293632046101 +0.5 21.31578947368421 0.006224821291274285 -0.006367629842832951 +0.8368421052631579 21.31578947368421 0.010524956121927553 -0.0065636378336877445 +1.5 21.31578947368421 0.01953772673322841 -0.007263814783676106 +1.6236842105263158 21.31578947368421 0.021334091794805922 -0.007444712650788156 +2.4105263157894736 21.31578947368421 0.03402551019114632 -0.009031936087471335 +2.5 21.31578947368421 0.035640008049888365 -0.009266400886054023 +3.1973684210526314 21.31578947368421 0.04985332443874577 -0.011562879746081858 +3.5 21.31578947368421 0.05712116073273154 -0.012864854560120239 +3.9842105263157896 21.31578947368421 0.07054911117439742 -0.015434690273531858 +4.771052631578947 21.31578947368421 0.09862524398279209 -0.021290407379179258 +5.557894736842106 21.31578947368421 0.13785128370930763 -0.030166061257600692 +6.344736842105263 21.31578947368421 0.19402432744416231 -0.04374119048488001 +7.131578947368421 21.31578947368421 0.27624109153171433 -0.06476909183705205 +7.918421052631579 21.31578947368421 0.3990354249729334 -0.09782469664862455 +8.705263157894738 21.31578947368421 0.5860279632796679 -0.15062234375465508 +9.492105263157894 21.31578947368421 0.8762550156672603 -0.2363711851522142 +10.278947368421052 21.31578947368421 1.3353137304860505 -0.3780475556007622 +11.06578947368421 21.31578947368421 2.075293872591344 -0.6162608839471911 +11.852631578947369 21.31578947368421 3.2909809017925142 -1.0239556404481231 +12.639473684210527 21.31578947368421 5.326648640361277 -1.734310909261443 +13.426315789473684 21.31578947368421 8.801239623440983 -2.9945000693017274 +14.213157894736842 21.31578947368421 14.846707414633615 -5.270882723432778 +15.0 21.31578947368421 25.5690680432627 -9.45802889239218 +0.05 22.105263157894736 0.0005867196865858582 -0.005932336167665591 +0.5 22.105263157894736 0.005898972934291435 -0.006029570933498933 +0.8368421052631579 22.105263157894736 0.009970509725079476 -0.006208722384653577 +1.5 22.105263157894736 0.01848616909911272 -0.00684786984590236 +1.6236842105263158 22.105263157894736 0.020179782466896598 -0.007012801344388906 +2.4105263157894736 22.105263157894736 0.0321052369360967 -0.008456721066685784 +2.5 22.105263157894736 0.03361715911970466 -0.008669559091371893 +3.1973684210526314 22.105263157894736 0.046879330711492494 -0.010748697109112786 +3.5 22.105263157894736 0.05362960284870931 -0.011923433794442741 +3.9842105263157896 22.105263157894736 0.06605211888636109 -0.014234754568790496 +4.771052631578947 22.105263157894736 0.09185209543076883 -0.019471120610297085 +5.557894736842106 22.105263157894736 0.12759171661463573 -0.027346555832593322 +6.344736842105263 22.105263157894736 0.17831650829648524 -0.03928873950631778 +7.131578947368421 22.105263157894736 0.25186589071982224 -0.05761513374023858 +7.918421052631579 22.105263157894736 0.3606367988903776 -0.08613631894876835 +8.705263157894738 22.105263157894736 0.524561845002215 -0.13120709709990938 +9.492105263157894 22.105263157894736 0.7762147889083888 -0.20358313273696105 +10.278947368421052 22.105263157894736 1.1696945004962458 -0.3217473843848041 +11.06578947368421 22.105263157894736 1.796316758975061 -0.5179565609534029 +11.852631578947369 22.105263157894736 2.8127508510242105 -0.8493934708454958 +12.639473684210527 22.105263157894736 4.492240472593227 -1.4190421220348133 +13.426315789473684 22.105263157894736 7.319280042462254 -2.4153456117624374 +14.213157894736842 22.105263157894736 12.167315667676572 -4.18867810951316 +15.0 22.105263157894736 20.637449795512417 -7.4010275868801765 +0.05 22.894736842105264 0.000557143440253961 -0.005631212855830784 +0.5 22.894736842105264 0.005600599851235396 -0.005720400801810851 +0.8368421052631579 22.894736842105264 0.00946309563764385 -0.005884657819471923 +1.5 22.894736842105264 0.01752562725755622 -0.006469970562048135 +1.6236842105263158 22.894736842105264 0.01912587273546146 -0.006620840123129374 +2.4105263157894736 22.894736842105264 0.03035853800976673 -0.007938912000239343 +2.5 22.894736842105264 0.031778090908176106 -0.008132809225563028 +3.1973684210526314 22.894736842105264 0.04418764815370574 -0.010022231429152135 +3.5 22.894736842105264 0.050476511662196494 -0.01108636562102093 +3.9842105263157896 22.894736842105264 0.062006787090658384 -0.0131738325978189 +4.771052631578947 22.894736842105264 0.08580216330955069 -0.01787753515096983 +5.557894736842106 22.894736842105264 0.11849997641595957 -0.02490039238556898 +6.344736842105263 22.894736842105264 0.1645170793279358 -0.035464001391087185 +7.131578947368421 22.894736842105264 0.23065210970475508 -0.05153295011411821 +7.918421052631579 22.894736842105264 0.32755209613017705 -0.07630578462734976 +8.705263157894738 22.894736842105264 0.4721632843643598 -0.11506122630622014 +9.492105263157894 22.894736842105264 0.6918871796355983 -0.17663625894259216 +10.278947368421052 22.894736842105264 1.031731167101415 -0.27604363392074205 +11.06578947368421 22.894736842105264 1.5667912434107047 -0.439172461146272 +11.852631578947369 22.894736842105264 2.4243632538491786 -0.7113500826038125 +12.639473684210527 22.894736842105264 3.8236892760975696 -1.1731615609919865 +13.426315789473684 22.894736842105264 6.148468976156693 -1.9701026025497925 +14.213157894736842 22.894736842105264 10.081128664394187 -3.3689730811766774 +15.0 22.894736842105264 16.85519170512988 -5.866696141651443 +0.05 23.68421052631579 0.0005299718327554232 -0.005354734020700309 +0.5 23.68421052631579 0.0053265639338897506 -0.005436777672057273 +0.8368421052631579 23.68421052631579 0.008997311242731443 -0.0055878181594848095 +1.5 23.68421052631579 0.016645446120436908 -0.006125436177184528 +1.6236842105263158 23.68421052631579 0.018160555197068808 -0.006263865958420084 +2.4105263157894736 23.68421052631579 0.028764266164834244 -0.0074709074750339876 +2.5 23.68421052631579 0.030100309137543857 -0.007648136649632569 +3.1973684210526314 23.68421052631579 0.041742343725210435 -0.009371127811238241 +3.5 23.68421052631579 0.04761800730574466 -0.010338614903725658 +3.9842105263157896 23.68421052631579 0.058352738194672994 -0.012231179088985796 +4.771052631578947 23.68421052631579 0.08037361271409893 -0.016474077495054692 +5.557894736842106 23.68421052631579 0.11040262765796938 -0.02276555834906203 +6.344736842105263 23.68421052631579 0.1523270304597564 -0.03215726717743982 +7.131578947368421 23.68421052631579 0.21207731929536527 -0.04632567586173363 +7.918421052631579 23.68421052631579 0.2988557610571085 -0.06797483032893649 +8.705263157894738 23.68421052631579 0.42716907605482163 -0.10152349802447233 +9.492105263157894 23.68421052631579 0.6202400780639438 -0.15429273685862005 +10.278947368421052 23.68421052631579 0.9158150136697413 -0.23858582436051595 +11.06578947368421 23.68421052631579 1.3761891935337844 -0.3753808512135844 +11.852631578947369 23.68421052631579 2.1057619352441033 -0.6009770882959362 +12.639473684210527 23.68421052631579 3.282217074595082 -0.9791250963428167 +13.426315789473684 23.68421052631579 5.21270615406652 -1.6234768473950705 +14.213157894736842 23.68421052631579 8.436557297294673 -2.739724648478601 +15.0 23.68421052631579 13.915786681141169 -4.705822265551006 +0.05 24.473684210526315 0.0005049397232781127 -0.00510016509978279 +0.5 24.473684210526315 0.00507417341003577 -0.005175840862376849 +0.8368421052631579 24.473684210526315 0.008568526537468518 -0.005315107265863164 +1.5 24.473684210526315 0.015836521856225777 -0.005810298774539173 +1.6236842105263158 24.473684210526315 0.017273748823512042 -0.0059376781400214415 +2.4105263157894736 24.473684210526315 0.027304468629506657 -0.007046339407308028 +2.5 24.473684210526315 0.028564731808673294 -0.007208835032531391 +3.1973684210526314 24.473684210526315 0.03951313086664693 -0.00878514248002012 +3.5 24.473684210526315 0.04501721872459386 -0.00966776473525378 +3.9842105263157896 24.473684210526315 0.05503950233313941 -0.011389769128947604 +4.771052631578947 24.473684210526315 0.07548209114311887 -0.015231853802120117 +5.557894736842106 24.473684210526315 0.10315740291198558 -0.020892214429436617 +6.344736842105263 24.473684210526315 0.1415035974445027 -0.029281318674806923 +7.131578947368421 24.473684210526315 0.19572190953549695 -0.041838556351412204 +7.918421052631579 24.473684210526315 0.27381270423776727 -0.06086509168525842 +8.705263157894738 24.473684210526315 0.38827362487854855 -0.0900863096460003 +9.492105263157894 24.473684210526315 0.5589220029122435 -0.1356141262346217 +10.278947368421052 24.473684210526315 0.817650626341633 -0.20761481523151004 +11.06578947368421 24.473684210526315 1.2165517546982763 -0.32323791301106003 +11.852631578947369 24.473684210526315 1.8419873011967398 -0.511828851385965 +12.639473684210527 24.473684210526315 2.839297516669794 -0.8243319485409416 +13.426315789473684 24.473684210526315 4.456802119303987 -1.3504845087311348 +14.213157894736842 24.473684210526315 7.12526650781937 -2.250690803733763 +15.0 24.473684210526315 11.603461243783999 -3.8159278028695343 +0.05 25.263157894736842 0.0004818184177989134 -0.004865149707398226 +0.5 25.263157894736842 0.004841107619747338 -0.004935128761150816 +0.8368421052631579 25.263157894736842 0.008172753241420729 -0.005063867233624949 +1.5 25.263157894736842 0.015091034601470927 -0.005521175301537476 +1.6236842105263158 25.263157894736842 0.01645680030549193 -0.0056387003727836095 +2.4105263157894736 25.263157894736842 0.025963817505418475 -0.006659838988631236 +2.5 25.263157894736842 0.02715507841045749 -0.0068092567600956545 +3.1973684210526314 25.263157894736842 0.037474325889854755 -0.008255723860435412 +3.5 25.263157894736842 0.04264296911191843 -0.009063488338739626 +3.9842105263157896 25.263157894736842 0.05202461379728243 -0.01063552846020441 +4.771052631578947 25.263157894736842 0.07105723280167786 -0.01412721166052784 +5.557894736842106 25.263157894736842 0.0966467152364772 -0.019239969812378317 +6.344736842105263 25.263157894736842 0.13184807764871725 -0.026766199290192288 +7.131578947368421 25.263157894736842 0.18124592009718044 -0.03794878427897102 +7.918421052631579 25.263157894736842 0.2518336802507406 -0.05475809953446665 +8.705263157894738 25.263157894736842 0.354441898130806 -0.08035581365519791 +9.492105263157894 25.263157894736842 0.5060900805442556 -0.11988086479017528 +10.278947368421052 25.263157894736842 0.7339114977563378 -0.18179816210102834 +11.06578947368421 25.263157894736842 1.0817907673610774 -0.2802426796359701 +11.852631578947369 25.263157894736842 1.6217405753617944 -0.4391467752970427 +12.639473684210527 25.263157894736842 2.473665864696513 -0.6996032101299975 +13.426315789473684 25.263157894736842 3.8401686229407757 -1.1331754412706125 +14.213157894736842 25.263157894736842 6.068683127515592 -1.8662815357736091 +15.0 25.263157894736842 9.763939611370313 -3.1254669043974754 +0.05 26.05263157894737 0.00046040972617741345 -0.004647647612019377 +0.5 26.05263157894737 0.004625356331098745 -0.004712512833047592 +0.8368421052631579 26.05263157894737 0.0078065393299885785 -0.0048318049250197555 +1.5 26.05263157894737 0.014402233930972182 -0.005255165825943982 +1.6236842105263158 26.05263157894737 0.01570224456915966 -0.00536387133278571 +2.4105263157894736 26.05263157894737 0.024729156244526686 -0.0063068530290280485 +2.5 26.05263157894737 0.025857384045868976 -0.006444617335625967 +3.1973684210526314 26.05263157894737 0.03560402415317681 -0.0077756874710502265 +3.5 26.05263157894737 0.040468741905086855 -0.008517140944573504 +3.9842105263157896 26.05263157894737 0.04927212265638585 -0.009956742610126855 +4.771052631578947 26.05263157894737 0.06703995313686871 -0.01314064805473008 +5.557894736842106 26.05263157894737 0.09077269067836444 -0.01777583632445143 +6.344736842105263 26.05263157894737 0.12319660383660302 -0.02455533279336728 +7.131578947368421 26.05263157894737 0.16837169643357092 -0.034558045173702846 +7.918421052631579 26.05263157894737 0.23244227986441848 -0.049480785354263306 +8.705263157894738 26.05263157894737 0.32484581884677055 -0.0720233843679553 +9.492105263157894 26.05263157894737 0.4602859225644195 -0.10653539804305255 +10.278947368421052 26.05263157894737 0.6619947290891236 -0.1601153455572948 +11.06578947368421 26.05263157894737 0.9671977862656714 -0.2445024940997507 +11.852631578947369 26.05263157894737 1.4363884233509807 -0.37937393973396333 +12.639473684210527 26.05263157894737 2.169274988952974 -0.5981646759643277 +13.426315789473684 26.05263157894737 3.3325668883956956 -0.9584743958141603 +14.213157894736842 26.05263157894737 5.20903704010941 -1.5609186816459657 +15.0 26.05263157894737 8.285321107509558 -2.583725118149001 +0.05 26.842105263157894 0.0004405411316696773 -0.004445884388038304 +0.5 26.842105263157894 0.004425170445994223 -0.004506144145771282 +0.8368421052631579 26.842105263157894 0.007466883459389499 -0.004616932621917936 +1.5 26.842105263157894 0.013764265428216182 -0.005009772020886896 +1.6236842105263158 26.842105263157894 0.015003611343125325 -0.005110557131101866 +2.4105263157894736 26.842105263157894 0.02358913577632235 -0.005983498644808922 +2.5 26.842105263157894 0.02465960993983806 -0.006110840806678877 +3.1973684210526314 26.842105263157894 0.033883444442312255 -0.007338961394879797 +3.5 26.842105263157894 0.03847186069371162 -0.008021441542365941 +3.9842105263157896 26.842105263157894 0.04675142096971592 -0.0093435990725406 +4.771052631578947 26.842105263157894 0.06338033582132331 -0.012255972040035598 +5.557894736842106 26.842105263157894 0.0854533278634549 -0.016472676425038087 +6.344736842105263 26.842105263157894 0.11541308638487785 -0.022602611040109338 +7.131578947368421 26.842105263157894 0.1568707666981205 -0.03158698664992823 +7.918421052631579 26.842105263157894 0.21525023558472037 -0.04489485130133864 +8.705263157894738 26.842105263157894 0.2988172301557275 -0.06484494060787428 +9.492105263157894 26.842105263157894 0.42034493443311416 -0.09514146946494005 +10.278947368421052 26.842105263157894 0.5998439974228416 -0.14177664241656654 +11.06578947368421 26.842105263157894 0.8690941997344016 -0.21456934775206252 +11.852631578947369 26.842105263157894 1.2792628007353561 -0.3298208771737645 +12.639473684210527 26.842105263157894 1.9138766618835281 -0.5149558644816598 +13.426315789473684 26.842105263157894 2.9111966400723235 -0.8167347636656341 +14.213157894736842 26.842105263157894 4.503312743810295 -1.3159712134436965 +15.0 26.842105263157894 7.085353384385173 -2.1542442147757663 +0.05 27.63157894736842 0.0004220618392912132 -0.004258310277582318 +0.5 27.63157894736842 0.0042390216983100755 -0.004314409754387353 +0.8368421052631579 27.63157894736842 0.007151165076245121 -0.0044175197573793745 +1.5 27.63157894736842 0.013172029533082939 -0.004782831384863694 +1.6236842105263158 27.63157894736842 0.0143552678557787 -0.004876480788358035 +2.4105263157894736 27.63157894736842 0.022533920467200555 -0.005686447389829766 +2.5 27.63157894736842 0.023551328953419805 -0.005804436655140667 +3.1973684210526314 27.63157894736842 0.03229640324587868 -0.006940385362833209 +3.5 27.63157894736842 0.036632833728711206 -0.007570222626229005 +3.9842105263157896 27.63157894736842 0.04443630963696809 -0.008787829492705785 +4.771052631578947 27.63157894736842 0.06003596837008817 -0.011459656820519432 +5.557894736842106 27.63157894736842 0.08061950216034035 -0.015308014429045202 +6.344736842105263 27.63157894736842 0.10838376283063957 -0.020870186468378762 +7.131578947368421 27.63157894736842 0.14655381242644808 -0.02897107091094364 +7.918421052631579 27.63157894736842 0.1999387509553368 -0.040888886469492514 +8.705263157894738 27.63157894736842 0.2758127254208657 -0.058625783685873886 +9.492105263157894 27.63157894736842 0.3853292730074275 -0.0853546160667714 +10.278947368421052 27.63157894736842 0.5458202146584081 -0.12616502916337002 +11.06578947368421 27.63157894736842 0.7845786710105528 -0.1893241162527606 +11.852631578947369 27.63157894736842 1.1451617762240234 -0.2884321303417995 +12.639473684210527 27.63157894736842 1.6980228556540804 -0.4461535982995208 +13.426315789473684 27.63157894736842 2.558673425773122 -0.7007545232680092 +14.213157894736842 27.63157894736842 3.9191022787814487 -1.1176980178635534 +15.0 27.63157894736842 6.102824606719452 -1.8104692152572686 +0.05 28.42105263157895 0.000404839523532659 -0.004083566371478256 +0.5 28.42105263157895 0.004065569503060972 -0.004135896904341935 +0.8368421052631579 28.42105263157895 0.006857086985571054 -0.00423205340737186 +1.5 28.42105263157895 0.012621065930155209 -0.00457246379980025 +1.6236842105263158 28.42105263157895 0.013752290096281768 -0.004659665025305459 +2.4105263157894736 28.42105263157895 0.02155494891113678 -0.005412832190166483 +2.5 28.42105263157895 0.02252346994033396 -0.005522401026028992 +3.1973684210526314 28.42105263157895 0.030828890188838533 -0.006575550932628778 +3.5 28.42105263157895 0.03493482618192391 -0.007158231855312587 +3.9842105263157896 28.42105263157895 0.04230425096103803 -0.008282427733120993 +4.771052631578947 28.42105263157895 0.05697062067905155 -0.010740334081235614 +5.557894736842106 28.42105263157895 0.07621260897651837 -0.014263117711746839 +6.344736842105263 28.42105263157895 0.10201295123232185 -0.019326782873391176 +7.131578947368421 28.42105263157895 0.13726293167560832 -0.026657433385157685 +7.918421052631579 28.42105263157895 0.1862442446287707 -0.03737245830443003 +8.705263157894738 28.42105263157895 0.2553870770606543 -0.05320935659723033 +9.492105263157894 28.42105263157895 0.3544777439348512 -0.076900536929323 +10.278947368421052 28.42105263157895 0.4986059440330372 -0.11279407266626328 +11.06578947368421 28.42105263157895 0.7113426320384496 -0.16789363186831993 +11.852631578947369 28.42105263157895 1.029989141548845 -0.2536210273341026 +12.639473684210527 28.42105263157895 1.5143533935917304 -0.3888389413367883 +13.426315789473684 28.42105263157895 2.261603696643512 -0.6050969274362425 +14.213157894736842 28.42105263157895 3.4317206403227756 -0.9558459846032974 +15.0 28.42105263157895 5.2916534808777795 -1.5328209361853342 +0.05 29.210526315789473 0.00038875763638379813 -0.003920456646752127 +0.5 29.210526315789473 0.003903633531284545 -0.003969363479389354 +0.8368421052631579 29.210526315789473 0.006582627883007124 -0.004059205758624997 +1.5 29.210526315789473 0.012107458290947476 -0.004377027836732044 +1.6236842105263158 29.210526315789473 0.013190356819185185 -0.0044583855567460855 +2.4105263157894736 29.210526315789473 0.020644738097116956 -0.005160172082431776 +2.5 29.210526315789473 0.021568108614639432 -0.005262136936908517 +3.1973684210526314 29.210526315789473 0.02946872286887685 -0.0062406734365785185 +3.5 29.210526315789473 0.033363232357783386 -0.006780973697333952 +3.9842105263157896 29.210526315789473 0.04033576573122622 -0.007821426031215313 +4.771052631578947 29.210526315789473 0.05415318798796973 -0.01008839616580041 +5.557894736842106 29.210526315789473 0.07218269519481975 -0.013322280562547235 +6.344736842105263 29.210526315789473 0.0962197138960395 -0.01794639125365973 +7.131578947368421 29.210526315789473 0.12886561820589298 -0.024602481111318065 +7.918421052631579 29.210526315789473 0.17394736465796767 -0.03427164070869181 +8.705263157894738 29.210526315789473 0.2371729664761222 -0.04846882415462489 +9.492105263157894 29.210526315789473 0.327167975853038 -0.0695590618700252 +10.278947368421052 29.210526315789473 0.4571340054055413 -0.10127706185365097 +11.06578947368421 29.210526315789473 0.6475339751500271 -0.14959057127597045 +11.852631578947369 29.210526315789473 0.9304911509548666 -0.2241512840040975 +12.639473684210527 29.210526315789473 1.3570815736612205 -0.34076132707600615 +13.426315789473684 29.210526315789473 2.0095680200281625 -0.525615189775072 +14.213157894736842 29.210526315789473 3.0221722493005516 -0.8226811227517237 +15.0 29.210526315789473 4.616774374331995 -1.306698140329424 +0.05 30.0 0.0003737131671466522 -0.0037679247207267562 +0.5 30.0 0.0037521708989962576 -0.0038137134702657042 +0.8368421052631579 30.0 0.006326002908341638 -0.0038978071693535816 +1.5 30.0 0.011627755345366433 -0.004195084817292275 +1.6236842105263158 30.0 0.012665661780661886 -0.004271132729716904 +2.4105263157894736 30.0 0.019796722134350356 -0.00492631095779411 +2.5 30.0 0.02067829544396448 -0.005021389402516511 +3.1973684210526314 30.0 0.028205264483405396 -0.005932488684633226 +3.5 30.0 0.03190532670146649 -0.006434582116220524 +3.9842105263157896 30.0 0.03851394365342141 -0.007399716003592714 +4.771052631578947 30.0 0.05155683942234392 -0.009495680733419347 +5.557894736842106 30.0 0.06848696666447286 -0.012472261552190993 +6.344736842105263 30.0 0.09093521582720972 -0.016707254526885904 +7.131578947368421 30.0 0.12125003751219945 -0.022770040532551196 +7.918421052631579 30.0 0.16286444948763612 -0.03152559825890367 +8.705263157894738 30.0 0.2208653799515671 -0.04430070504167044 +9.492105263157894 30.0 0.302887589843305 -0.06315214978669359 +10.278947368421052 30.0 0.42053361230259045 -0.0913041389028175 +11.06578947368421 30.0 0.5916552913707827 -0.13386939663370565 +11.852631578947369 30.0 0.8440620905688163 -0.19905117429393596 +12.639473684210527 30.0 1.2216184548803695 -0.3001694846803214 +13.426315789473684 30.0 1.7943867630354546 -0.4591156856527174 +14.213157894736842 30.0 2.675699299227 -0.712310080790577 +15.0 30.0 4.051234553636601 -1.121095278894567 diff --git a/test/enzyme_test.jl b/test/enzyme_test.jl index 63f3abb..cdb55a4 100644 --- a/test/enzyme_test.jl +++ b/test/enzyme_test.jl @@ -11,6 +11,6 @@ for line in eachline("data/besselk/enzyme/besselkx_levin_enzyme_tests.csv") (v, x, dv, dx) = parse.(Float64, split(line)) test_dv = dbesselkx_dv(v, x) test_dx = dbesselkx_dx(v, x) - @test isapprox(dv, test_dv, rtol=1e-4) - @test isapprox(dx, test_dx, rtol=1e-4) + @test isapprox(dv, test_dv, rtol=5e-14) + @test isapprox(dx, test_dx, rtol=5e-14) end From 9c6d06313157e323fdaf4402f60986a74fe37e4c Mon Sep 17 00:00:00 2001 From: Chris Geoga Date: Wed, 10 May 2023 08:39:09 -0400 Subject: [PATCH 10/17] Update test README that shows how I generated the file --- test/data/besselk/enzyme/README.md | 9 +++++---- 1 file changed, 5 insertions(+), 4 deletions(-) diff --git a/test/data/besselk/enzyme/README.md b/test/data/besselk/enzyme/README.md index 33bcb9f..4732d58 100644 --- a/test/data/besselk/enzyme/README.md +++ b/test/data/besselk/enzyme/README.md @@ -1,13 +1,14 @@ ```julia -using FiniteDifferences, ArbNumerics, DelimitedFiles +using FiniteDifferences, ArbNumerics, DelimitedFiles, SpecialFunctions if !(@isdefined vx) - const FD = central_fdm(10,1) vgrid = sort(vcat(range(0.05, 15.0, length=20), [0.5, 1.5, 2.5, 3.5])) xgrid = range(15.0, 30.0, length=20) const vx = vec(collect(Iterators.product(vgrid, xgrid))) end +simplefd(f,x,h=ArbReal(1e-40)) = (f(x+h)-f(x))/h + ArbNumerics.setprecision(ArbReal, digits=100) function arb_besselkx(v,x) @@ -17,8 +18,8 @@ end ref_values = map(vx) do vxj (v,x) = vxj - dx = FD(_x->arb_besselkx(v, _x), x) - dv = FD(_v->arb_besselkx(_v, x), v) + dx = simplefd(_x->arb_besselkx(v, _x), x) + dv = simplefd(_v->arb_besselkx(_v, x), v) Float64.((dv, dx)) end From 6f3445c667dabc126f3808deef8bbdcb5e97ffda Mon Sep 17 00:00:00 2001 From: Chris Geoga Date: Wed, 10 May 2023 13:37:11 -0400 Subject: [PATCH 11/17] @heltonmc's nice convergence check that fixes the scalar evaluation and simplifies the manual rule. --- ext/BesselsEnzymeCoreExt.jl | 15 +++++++-------- src/Math/Math.jl | 2 +- 2 files changed, 8 insertions(+), 9 deletions(-) diff --git a/ext/BesselsEnzymeCoreExt.jl b/ext/BesselsEnzymeCoreExt.jl index d21493c..15f43a0 100644 --- a/ext/BesselsEnzymeCoreExt.jl +++ b/ext/BesselsEnzymeCoreExt.jl @@ -6,10 +6,9 @@ module BesselsEnzymeCoreExt using EnzymeCore.EnzymeRules using Bessels.Math - # A manual method that takes an NTuple of partial sum terms and checks if it is - # exactly converged before computing the Levin sequence transformation. If it is - # exactly converged, the sequence transformation will create a divide-by-zero - # problem. See + # A manual method that separately transforms the `val` and `dval`, because + # sometimes the `val` can converge while the `dval` hasn't, so just using an + # early return or something can give incorrect derivatives in edge cases. # # https://github.com/JuliaMath/Bessels.jl/issues/96 # @@ -30,10 +29,10 @@ module BesselsEnzymeCoreExt ::Type{<:Duplicated}, s::Duplicated, w::Duplicated) - (sv, dv, N) = (s.val, s.dval, length(s.val)) - ls = (sv[N-1] == sv[N]) ? sv[N] : levin_transform(sv, w.val) - dls = (dv[N-1] == dv[N]) ? dv[N] : levin_transform(dv, w.dval) - Duplicated(ls, dls) + (sv, dv, N) = (s.val, s.dval, length(s.val)) + ls = levin_transform(sv, w.val) + dls = levin_transform(dv, w.dval) + Duplicated(ls, dls) end end diff --git a/src/Math/Math.jl b/src/Math/Math.jl index 386ca64..73f608a 100644 --- a/src/Math/Math.jl +++ b/src/Math/Math.jl @@ -136,7 +136,7 @@ end len = N - 1 :( begin - @nexprs $N i -> a_{i} = Vec{2, T}((s[i] / w[i], 1 / w[i])) + @nexprs $N i -> a_{i} = iszero(w[i]) ? (return s[i]) : Vec{2, T}((s[i] / w[i], 1 / w[i])) @nexprs $len k -> (@nexprs ($len-k) i -> a_{i} = fmadd(a_{i}, levin_scale(one(T), i, k-1), a_{i+1})) return (a_1[1] / a_1[2]) end From 3c303947ff9bbfdc9e1385e86b45d166fe6a37bb Mon Sep 17 00:00:00 2001 From: Chris Geoga Date: Tue, 16 May 2023 18:20:11 -0400 Subject: [PATCH 12/17] Power series Temme-like method for near-integer/integer orders that gives valid derivatives. Tests TODO/to debug still. --- src/BesselFunctions/besselk.jl | 114 ++++++++++++++++++++++++--------- src/GammaFunctions/gamma.jl | 2 + src/Math/Math.jl | 3 + 3 files changed, 88 insertions(+), 31 deletions(-) diff --git a/src/BesselFunctions/besselk.jl b/src/BesselFunctions/besselk.jl index 7407464..abbb776 100644 --- a/src/BesselFunctions/besselk.jl +++ b/src/BesselFunctions/besselk.jl @@ -499,38 +499,25 @@ besselk_power_series(v, x::Float32) = Float32(besselk_power_series(v, Float64(x) besselk_power_series(v, x::ComplexF32) = ComplexF32(besselk_power_series(v, ComplexF64(x))) function besselk_power_series(v, x::ComplexOrReal{T}) where T - MaxIter = 1000 - S = eltype(x) - v, x = S(v), S(x) - - z = x / 2 - zz = z * z - logz = log(z) - xd2_v = exp(v*logz) - xd2_nv = inv(xd2_v) - - # use the reflection identify to calculate gamma(-v) - # use relation gamma(v)*v = gamma(v+1) to avoid two gamma calls - gam_v = gamma(v) - gam_nv = π / (sinpi(-abs(v)) * gam_v * v) - gam_1mv = -gam_nv * v - gam_1mnv = gam_v * v - - _t1 = gam_v * xd2_nv * gam_1mv - _t2 = gam_nv * xd2_v * gam_1mnv - (xd2_pow, fact_k, out) = (one(S), one(S), zero(S)) - for k in 0:MaxIter - t1 = xd2_pow * T(0.5) - tmp = muladd(_t1, gam_1mnv, _t2 * gam_1mv) - tmp *= inv(gam_1mv * gam_1mnv * fact_k) - term = t1 * tmp - out += term - abs(term / out) < eps(T) && break - (gam_1mnv, gam_1mv) = (gam_1mnv*(one(S) + v + k), gam_1mv*(one(S) - v + k)) - xd2_pow *= zz - fact_k *= k + one(S) + Math.isnearint(v) && return besselk_power_series_int(v, x) + MaxIter = 5000 + gam = gamma(v) + ngam = π / (sinpi(-abs(v)) * gam * v) + + s1, s2 = zero(T), zero(T) + t1, t2 = one(T), one(T) + + for k in 1:MaxIter + s1 += t1 + s2 += t2 + t1 *= x^2 / (4k * (k - v)) + t2 *= x^2 / (4k * (k + v)) + abs(t1) < eps(T) && break end - return out + + xpv = (x/2)^v + s = gam * s1 + xpv^2 * ngam * s2 + return s / (2*xpv) end besselk_power_series_cutoff(nu, x::Float64) = x < 2.0 || nu > 1.6x - 1.0 besselk_power_series_cutoff(nu, x::Float32) = x < 10.0f0 || nu > 1.65f0*x - 8.0f0 @@ -615,3 +602,68 @@ end end ) end + +# This is an expansion of the function +# +# f_0(v, x) = (x^v)*gamma(-v) + (x^(-v))*gamma(v) +# = (x^v)*(gamma(-v) + (x^(-2*v))*gamma(v)) +# +# around v ∼ 0. As you can see by plugging that second form into Wolfram Alpha +# and getting an expansion back, this is actually a bivariate polynomial in +# (v^2, log(x)). So that's how this is structured. +@inline function f0_local_expansion_v0(v, x) + lx = log(x) + c0 = evalpoly(lx, (-1.1544313298030657, -2.0)) + c2 = evalpoly(lx, ( 1.4336878944573288, -1.978111990655945, -0.5772156649015329, -0.3333333333333333)) + c4 = evalpoly(lx, (-0.6290784463642211, -1.4584260788225176, -0.23263776388631713, -0.32968533177599085, -0.048101305408461074, -0.016666666666666666)) + evalpoly(v*v, (c0,c2,c4))/2 +end + +# This function assumes |v| < 1e-6 or 1e-7! +# +# TODO (cg 2023/05/16 18:07): lots of micro-optimizations. +function besselk_power_series_temme_basal(v::V, x::Float64) where{V} + max_iter = 50 + T = promote_type(V,Float64) + z = x/2 + zz = z*z + fk = f0_local_expansion_v0(v, x/2) + zv = z^v + znv = inv(zv) + gam_1pv = GammaFunctions.gamma_near_1(1+v) + gam_1nv = GammaFunctions.gamma_near_1(1-v) + (pk, qk, _ck, factk, vv) = (znv*gam_1pv/2, zv*gam_1nv/2, one(T), one(T), v*v) + (out_v, out_vp1) = (zero(T), zero(T)) + for k in 1:max_iter + # add to the series: + ck = _ck/factk + term_v = ck*fk + term_vp1 = ck*(pk - (k-1)*fk) + out_v += term_v + out_vp1 += term_vp1 + # check for convergence: + ((abs(term_v) < eps(T)) && (abs(term_vp1) < eps(T))) && break + # otherwise, increment new quantities: + fk = (k*fk + pk + qk)/(k^2 - vv) + pk /= (k-v) + qk /= (k+v) + _ck *= zz + factk *= k + end + (out_v, out_vp1/z) +end + +function besselk_power_series_int(v, x::Float64) + v < zero(v) && return besselk_power_series_int(-v, x) + flv = Int(floor(v)) + _v = v - flv + (kv, kvp1) = besselk_power_series_temme_basal(_v, x) + abs(v) < 1/2 && return kv + twodx = 2/x + for _ in 1:(flv-1) + _v += 1 + (kv, kvp1) = (kvp1, muladd(twodx*_v, kvp1, kv)) + end + kvp1 +end + diff --git a/src/GammaFunctions/gamma.jl b/src/GammaFunctions/gamma.jl index 6350cf2..9a9df75 100644 --- a/src/GammaFunctions/gamma.jl +++ b/src/GammaFunctions/gamma.jl @@ -113,3 +113,5 @@ function gamma(n::Integer) n > 20 && return gamma(float(n)) @inbounds return Float64(factorial(n-1)) end + +gamma_near_1(x) = evalpoly(x-one(x), (1.0, -0.5772156649015329, 0.9890559953279725, -0.23263776388631713)) \ No newline at end of file diff --git a/src/Math/Math.jl b/src/Math/Math.jl index 73f608a..b7793aa 100644 --- a/src/Math/Math.jl +++ b/src/Math/Math.jl @@ -154,4 +154,7 @@ end ) end +# TODO (cg 2023/05/16 18:09): dispute this cutoff. +isnearint(x) = abs(x-round(x)) < 1e-7 + end From 0e687a1adf396090b9eea3a677e938ac3adb1108 Mon Sep 17 00:00:00 2001 From: Chris Geoga Date: Wed, 17 May 2023 16:38:27 -0400 Subject: [PATCH 13/17] Checkpoint commit with small tweaks. Still chasing down a small issue that makes the AD just a bit less accurate than it should be. --- ext/BesselsEnzymeCoreExt.jl | 13 ++++++++++++ src/BesselFunctions/besselk.jl | 13 ++++++------ src/Math/Math.jl | 2 +- test/besselk_enzyme_test.jl | 38 ++++++++++++++++++++++++++++++++++ test/runtests.jl | 2 +- 5 files changed, 60 insertions(+), 8 deletions(-) create mode 100644 test/besselk_enzyme_test.jl diff --git a/ext/BesselsEnzymeCoreExt.jl b/ext/BesselsEnzymeCoreExt.jl index 15f43a0..74f87b1 100644 --- a/ext/BesselsEnzymeCoreExt.jl +++ b/ext/BesselsEnzymeCoreExt.jl @@ -35,4 +35,17 @@ module BesselsEnzymeCoreExt Duplicated(ls, dls) end + # This is fixing a straight bug in Enzyme. + function EnzymeRules.forward(func::Const{typeof(sinpi)}, + ::Type{<:Duplicated}, + x::Duplicated) + Duplicated(sinpi(x.val), pi*cospi(x.val)) + end + + function EnzymeRules.forward(func::Const{typeof(sinpi)}, + ::Type{<:Const}, + x::Const) + sinpi(x.val) + end + end diff --git a/src/BesselFunctions/besselk.jl b/src/BesselFunctions/besselk.jl index abbb776..3d2ade7 100644 --- a/src/BesselFunctions/besselk.jl +++ b/src/BesselFunctions/besselk.jl @@ -654,16 +654,17 @@ function besselk_power_series_temme_basal(v::V, x::Float64) where{V} end function besselk_power_series_int(v, x::Float64) - v < zero(v) && return besselk_power_series_int(-v, x) - flv = Int(floor(v)) - _v = v - flv + v = abs(v) + (_v, flv) = modf(v) + if _v > 1/2 + (_v, flv) = (_v-one(_v), flv+1) + end (kv, kvp1) = besselk_power_series_temme_basal(_v, x) - abs(v) < 1/2 && return kv twodx = 2/x - for _ in 1:(flv-1) + for _ in 1:flv _v += 1 (kv, kvp1) = (kvp1, muladd(twodx*_v, kvp1, kv)) end - kvp1 + kv end diff --git a/src/Math/Math.jl b/src/Math/Math.jl index b7793aa..b07cf09 100644 --- a/src/Math/Math.jl +++ b/src/Math/Math.jl @@ -155,6 +155,6 @@ end end # TODO (cg 2023/05/16 18:09): dispute this cutoff. -isnearint(x) = abs(x-round(x)) < 1e-7 +isnearint(x) = abs(x-round(x)) < 1e-5 end diff --git a/test/besselk_enzyme_test.jl b/test/besselk_enzyme_test.jl new file mode 100644 index 0000000..77f95c2 --- /dev/null +++ b/test/besselk_enzyme_test.jl @@ -0,0 +1,38 @@ +using EnzymeCore, Enzyme +import Bessels.BesselFunctions: besselkx_levin +import Bessels.BesselFunctions: besselk_power_series + +dbesselkx_dv(v, x) = autodiff(Forward, _v->besselkx_levin(_v, x, Val(30)), + Duplicated, Duplicated(v, 1.0))[2] + +dbesselkx_dx(v, x) = autodiff(Forward, _x->besselkx_levin(v, _x, Val(30)), + Duplicated, Duplicated(x, 1.0))[2] + +#= +dbesselk_ps_dv(v, x) = autodiff(Forward, _v->besselk_power_series(_v, x), + Duplicated, Duplicated(v, 1.0))[2] + +dbesselk_ps_dx(v, x) = autodiff(Forward, _x->besselk_power_series(v, _x), + Duplicated, Duplicated(x, 1.0))[2] +=# + + +for line in eachline("data/besselk/enzyme/besselkx_levin_enzyme_tests.csv") + (v, x, dv, dx) = parse.(Float64, split(line)) + test_dv = dbesselkx_dv(v, x) + test_dx = dbesselkx_dx(v, x) + @test isapprox(dv, test_dv, rtol=5e-14) + @test isapprox(dx, test_dx, rtol=5e-14) +end + +#= +for line in eachline("data/besselk/enzyme/besselk_power_series_enzyme_tests.csv") + (v, x, dv, dx) = parse.(Float64, split(line)) + test_dv = dbesselk_ps_dv(v, x) + test_dx = dbesselk_ps_dx(v, x) + #@test isapprox(dv, test_dv, rtol=5e-14) + #@test isapprox(dx, test_dx, rtol=5e-14) + +end +=# + diff --git a/test/runtests.jl b/test/runtests.jl index db22996..0a8932f 100644 --- a/test/runtests.jl +++ b/test/runtests.jl @@ -10,4 +10,4 @@ import SpecialFunctions @time @testset "gamma" begin include("gamma_test.jl") end @time @testset "airy" begin include("airy_test.jl") end @time @testset "sphericalbessel" begin include("sphericalbessel_test.jl") end -@time @testset "enzyme autodiff" begin include("enzyme_test.jl") end +@time @testset "besselk enzyme autodiff" begin include("besselk_enzyme_test.jl") end From 232f566d8938c44d55f8cf03615615a7da4c7ec1 Mon Sep 17 00:00:00 2001 From: Chris Geoga Date: Thu, 18 May 2023 11:36:18 -0400 Subject: [PATCH 14/17] near-int autodiff working, just with slightly less precision with 1e-8 of v=0. But probably ready to merge. --- src/BesselFunctions/besselk.jl | 39 ++-- test/besselk_enzyme_test.jl | 30 +-- .../besselk_power_series_enzyme_tests.csv | 220 ++++++++++++++++++ 3 files changed, 253 insertions(+), 36 deletions(-) create mode 100644 test/data/besselk/enzyme/besselk_power_series_enzyme_tests.csv diff --git a/src/BesselFunctions/besselk.jl b/src/BesselFunctions/besselk.jl index 3d2ade7..52ab073 100644 --- a/src/BesselFunctions/besselk.jl +++ b/src/BesselFunctions/besselk.jl @@ -603,35 +603,32 @@ end ) end -# This is an expansion of the function -# -# f_0(v, x) = (x^v)*gamma(-v) + (x^(-v))*gamma(v) -# = (x^v)*(gamma(-v) + (x^(-2*v))*gamma(v)) -# -# around v ∼ 0. As you can see by plugging that second form into Wolfram Alpha -# and getting an expansion back, this is actually a bivariate polynomial in -# (v^2, log(x)). So that's how this is structured. +# This is a version of Temme's proposed f_0 (1975 JCP, see reference above) that +# swaps in a bunch of local expansions for functions that are well-behaved but +# whose standard forms can't be naively evaluated by a computer at the origin. @inline function f0_local_expansion_v0(v, x) - lx = log(x) - c0 = evalpoly(lx, (-1.1544313298030657, -2.0)) - c2 = evalpoly(lx, ( 1.4336878944573288, -1.978111990655945, -0.5772156649015329, -0.3333333333333333)) - c4 = evalpoly(lx, (-0.6290784463642211, -1.4584260788225176, -0.23263776388631713, -0.32968533177599085, -0.048101305408461074, -0.016666666666666666)) - evalpoly(v*v, (c0,c2,c4))/2 + l2dx = log(2/x) + mu = v*l2dx + vv = v*v + sp = evalpoly(vv, (1.0, 1.6449340668482264, 1.8940656589944918, 1.9711021825948702)) + g1 = evalpoly(vv, (-0.5772156649015329, 0.04200263503409518, 0.042197734555544306)) + g2 = evalpoly(vv, (1.0, -0.6558780715202539, 0.16653861138229145)) + sh = evalpoly(mu*mu, (1.0, 0.16666666666666666, 0.008333333333333333, 0.0001984126984126984, 2.7557319223985893e-6)) + sp*(g1*cosh(mu) + g2*sh*l2dx) end -# This function assumes |v| < 1e-6 or 1e-7! -# -# TODO (cg 2023/05/16 18:07): lots of micro-optimizations. -function besselk_power_series_temme_basal(v::V, x::Float64) where{V} +# This function assumes |v|<1e-5! +function besselk_power_series_temme_basal(v::V, x::X) where{V,X} max_iter = 50 - T = promote_type(V,Float64) + T = promote_type(V,X) z = x/2 zz = z*z - fk = f0_local_expansion_v0(v, x/2) + fk = f0_local_expansion_v0(v,x) zv = z^v znv = inv(zv) - gam_1pv = GammaFunctions.gamma_near_1(1+v) - gam_1nv = GammaFunctions.gamma_near_1(1-v) + gam_1_c = (1.0, -0.5772156649015329, 0.9890559953279725, -0.23263776388631713) + gam_1pv = evalpoly(v, gam_1_c) + gam_1nv = evalpoly(-v, gam_1_c) (pk, qk, _ck, factk, vv) = (znv*gam_1pv/2, zv*gam_1nv/2, one(T), one(T), v*v) (out_v, out_vp1) = (zero(T), zero(T)) for k in 1:max_iter diff --git a/test/besselk_enzyme_test.jl b/test/besselk_enzyme_test.jl index 77f95c2..89fdcfb 100644 --- a/test/besselk_enzyme_test.jl +++ b/test/besselk_enzyme_test.jl @@ -8,31 +8,31 @@ dbesselkx_dv(v, x) = autodiff(Forward, _v->besselkx_levin(_v, x, Val(30)), dbesselkx_dx(v, x) = autodiff(Forward, _x->besselkx_levin(v, _x, Val(30)), Duplicated, Duplicated(x, 1.0))[2] -#= dbesselk_ps_dv(v, x) = autodiff(Forward, _v->besselk_power_series(_v, x), Duplicated, Duplicated(v, 1.0))[2] dbesselk_ps_dx(v, x) = autodiff(Forward, _x->besselk_power_series(v, _x), Duplicated, Duplicated(x, 1.0))[2] -=# for line in eachline("data/besselk/enzyme/besselkx_levin_enzyme_tests.csv") - (v, x, dv, dx) = parse.(Float64, split(line)) - test_dv = dbesselkx_dv(v, x) - test_dx = dbesselkx_dx(v, x) - @test isapprox(dv, test_dv, rtol=5e-14) - @test isapprox(dx, test_dx, rtol=5e-14) + (v, x, dv, dx) = parse.(Float64, split(line)) + test_dv = dbesselkx_dv(v, x) + test_dx = dbesselkx_dx(v, x) + @test isapprox(dv, test_dv, rtol=5e-14) + @test isapprox(dx, test_dx, rtol=5e-14) end -#= for line in eachline("data/besselk/enzyme/besselk_power_series_enzyme_tests.csv") - (v, x, dv, dx) = parse.(Float64, split(line)) - test_dv = dbesselk_ps_dv(v, x) - test_dx = dbesselk_ps_dx(v, x) - #@test isapprox(dv, test_dv, rtol=5e-14) - #@test isapprox(dx, test_dx, rtol=5e-14) - + (v, x, dv, dx) = parse.(Float64, split(line)) + test_dv = dbesselk_ps_dv(v, x) + test_dx = dbesselk_ps_dx(v, x) + if abs(v) <= 1e-8 + @test isapprox(dv, test_dv, rtol=1e-7) + @test isapprox(dx, test_dx, rtol=1e-7) + else + @test isapprox(dv, test_dv, rtol=5e-14) + @test isapprox(dx, test_dx, rtol=5e-14) + end end -=# diff --git a/test/data/besselk/enzyme/besselk_power_series_enzyme_tests.csv b/test/data/besselk/enzyme/besselk_power_series_enzyme_tests.csv new file mode 100644 index 0000000..514a103 --- /dev/null +++ b/test/data/besselk/enzyme/besselk_power_series_enzyme_tests.csv @@ -0,0 +1,220 @@ +-1.0e-8 1.0e-5 -5.425186613586617e-6 -99999.99993935639 +1.0e-8 1.0e-5 5.425186613586617e-6 -99999.99993935639 +0.99999999 1.0e-5 1.162885561219195e6 -9.999998742678802e9 +1.0 1.0e-5 1.1628856980944362e6 -1.0000000005564426e10 +1.00000001 1.0e-5 1.1628858349696922e6 -1.0000001268450197e10 +1.99999999 1.0e-5 2.5257710758016803e11 -3.9999994748457575e15 +2.0 1.0e-5 2.5257713960675833e11 -3.999999999999999e15 +2.00000001 1.0e-5 2.5257717163335272e11 -4.0000005251543095e15 +2.99999999 1.0e-5 1.0503084202282163e17 -2.399999676897455e21 +3.0 1.0e-5 1.050308558437662e17 -2.399999999989999e21 +3.00000001 1.0e-5 1.0503086966471262e17 -2.4000003230825867e21 +-1.0e-8 0.07895684210526316 -9.842689968784812e-8 -12.540480739771382 +1.0e-8 0.07895684210526316 9.842689968784812e-8 -12.540480739771382 +0.99999999 0.07895684210526316 33.695420311958394 -161.48750833859165 +1.0 0.07895684210526316 33.695421405878655 -161.48751419443647 +1.00000001 0.07895684210526316 33.69542249979894 -161.48752005028146 +1.99999999 0.07895684210526316 1171.1689358673987 -8126.201626601186 +2.0 0.07895684210526316 1171.1689807454445 -8126.20196416698 +2.00000001 0.07895684210526316 1171.1690256234922 -8126.20230173279 +2.99999999 0.07895684210526316 67479.46254161856 -617361.0529533329 +3.0 0.07895684210526316 67479.46540958171 -617361.0806609665 +3.00000001 0.07895684210526316 67479.46827754499 -617361.1083686013 +-1.0e-8 0.1579036842105263 -5.0440863814507895e-8 -6.137827365309258 +1.0e-8 0.1579036842105263 5.0440863814507895e-8 -6.137827365309258 +0.99999999 0.1579036842105263 12.540534440143462 -40.85089943450172 +1.0 0.1579036842105263 12.54053477956886 -40.85090061739761 +1.00000001 0.1579036842105263 12.540535118994264 -40.85090180029353 +1.99999999 0.1579036842105263 236.5791734885244 -1015.8875645600523 +2.0 0.1579036842105263 236.57918101486533 -1015.8875996992039 +2.00000001 0.1579036842105263 236.57918854120652 -1015.8876348383567 +2.99999999 0.1579036842105263 7015.289705250685 -38564.64939708375 +3.0 0.1579036842105263 7015.289956176605 -38564.650860562164 +3.00000001 0.1579036842105263 7015.290207102535 -38564.65232404064 +-1.0e-8 0.23685052631578948 -3.1868782720708587e-8 -3.9762225988287523 +1.0e-8 0.23685052631578948 3.1868782720708587e-8 -3.9762225988287523 +0.99999999 0.23685052631578948 6.7226297328406535 -18.380156491460472 +1.0 0.23685052631578948 6.722629895414007 -18.380156943173738 +1.00000001 0.23685052631578948 6.722630057987364 -18.380157394887007 +1.99999999 0.23685052631578948 90.34265078371085 -300.9403452225231 +2.0 0.23685052631578948 90.34265332370957 -300.94035440323375 +2.00000001 0.23685052631578948 90.34265586370834 -300.94036358394476 +2.99999999 0.23685052631578948 1829.4194890184042 -7608.355353997639 +3.0 0.23685052631578948 1829.4195473345358 -7608.355611863256 +3.00000001 0.23685052631578948 1829.4196056506692 -7608.355869728881 +-1.0e-8 0.3157973684210526 -2.2114994884074176e-8 -2.8823485371413677 +1.0e-8 0.3157973684210526 2.2114994884074176e-8 -2.8823485371413677 +0.99999999 0.3157973684210526 4.197558567860136 -10.45278660893378 +1.0 0.3157973684210526 4.197558661407828 -10.452786833125234 +1.00000001 0.3157973684210526 4.197558754955521 -10.45278705731669 +1.99999999 0.3157973684210526 44.83829015394938 -126.88588429915002 +2.0 0.3157973684210526 44.838291300195735 -126.88588780083069 +2.00000001 0.3157973684210526 44.83829244644212 -126.88589130251145 +2.99999999 0.3157973684210526 696.1385689986124 -2402.9698935262886 +3.0 0.3157973684210526 696.1385892902923 -2402.969968050816 +3.00000001 0.3157973684210526 696.1386095819729 -2402.970042575346 +-1.0e-8 0.3947442105263158 -1.620723523542244e-8 -2.2193650236523337 +1.0e-8 0.3947442105263158 1.620723523542244e-8 -2.2193650236523337 +0.99999999 0.3947442105263158 2.852736516334018 -6.74838756294082 +1.0 0.3947442105263158 2.852736575983206 -6.748387691431661 +1.00000001 0.3947442105263158 2.8527366356323944 -6.748387819922504 +1.99999999 0.3947442105263158 25.698167798180446 -64.89627569640321 +2.0 0.3947442105263158 25.698168405676494 -64.89627734033134 +2.00000001 0.3947442105263158 25.698169013172556 -64.89627898425951 +2.99999999 0.3947442105263158 325.93288572340526 -981.9087892240608 +3.0 0.3947442105263158 325.9328945429319 -981.9088174832733 +3.00000001 0.3947442105263158 325.9329033624587 -981.9088457424866 +-1.0e-8 0.47369105263157896 -1.2314475707383088e-8 -1.774033510508262 +1.0e-8 0.47369105263157896 1.2314475707383088e-8 -1.774033510508262 +0.99999999 0.47369105263157896 2.046729861500021 -4.714645320359511 +1.0 0.47369105263157896 2.0467299021275345 -4.714645401018907 +1.00000001 0.47369105263157896 2.0467299427550483 -4.714645481678303 +1.99999999 0.47369105263157896 16.131878813969152 -37.49255892821751 +2.0 0.47369105263157896 16.13187917065201 -37.49255980839138 +2.00000001 0.47369105263157896 16.131879527334874 -37.49256068856526 +2.99999999 0.47369105263157896 173.9880415637526 -472.12327507882605 +3.0 0.47369105263157896 173.98804597855303 -472.1232878047734 +3.00000001 0.47369105263157896 173.98805039335357 -472.12330053072105 +-1.0e-8 0.5526378947368421 -9.601672928833976e-9 -1.4544519301218142 +1.0e-8 0.5526378947368421 9.601672928833976e-9 -1.4544519301218142 +0.99999999 0.5526378947368421 1.5249033918939678 -3.4745546815205444 +1.0 0.5526378947368421 1.5249034208702144 -3.4745547354320725 +1.00000001 0.5526378947368421 1.5249034498464606 -3.474554789343601 +1.99999999 0.5526378947368421 10.782305438114463 -23.553515241792073 +2.0 0.5526378947368421 10.78230566295405 -23.553515757748738 +2.00000001 0.5526378947368421 10.782305887793642 -23.55351627370541 +2.99999999 0.5526378947368421 101.66642817255257 -253.93184048861562 +3.0 0.5526378947368421 101.66643060934518 -253.93184694149514 +3.00000001 0.5526378947368421 101.66643304613785 -253.93185339437483 +-1.0e-8 0.6315847368421053 -7.633873918032207e-9 -1.2143621871586956 +1.0e-8 0.6315847368421053 7.633873918032207e-9 -1.2143621871586956 +0.99999999 0.6315847368421053 1.1681623925342115 -2.6605160177707723 +1.0 0.6315847368421053 1.1681624139161837 -2.6605160554937326 +1.00000001 0.6315847368421053 1.1681624352981559 -2.6605160932166934 +1.99999999 0.6315847368421053 7.5445919286506316 -15.727817597972695 +2.0 0.6315847368421053 7.544592077976403 -15.727817921131486 +2.00000001 0.6315847368421053 7.544592227302177 -15.727818244290283 +2.99999999 0.6315847368421053 63.463595195347054 -148.22796810980893 +3.0 0.6315847368421053 63.46359664027001 -148.22797167856376 +3.00000001 0.6315847368421053 63.463598085193006 -148.22797524731868 +-1.0e-8 0.7105315789473684 -6.162808887534923e-9 -1.0278688932014883 +1.0e-8 0.7105315789473684 6.162808887534923e-9 -1.0278688932014883 +0.99999999 0.7105315789473684 0.9142131788593283 -2.0961968960798703 +1.0 0.7105315789473684 0.9142131950562388 -2.0961969234126747 +1.00000001 0.7105315789473684 0.9142132112531494 -2.0961969507454796 +1.99999999 0.7105315789473684 5.46656083735013 -11.000167187018414 +2.0 0.7105315789473684 5.466560940570338 -11.000167399894464 +2.00000001 0.7105315789473684 5.466561043790549 -11.000167612770516 +2.99999999 0.7105315789473684 41.66099628077927 -92.0925700744845 +3.0 0.7105315789473684 41.66099718580828 -92.09257218332277 +3.00000001 0.7105315789473684 41.660998090837296 -92.09257429216109 +-1.0e-8 0.7894784210526316 -5.036921200254442e-9 -0.8793030336266766 +1.0e-8 0.7894784210526316 5.036921200254442e-9 -0.8793030336266766 +0.99999999 0.7894784210526316 0.7277034283442436 -1.6882833030691788 +1.0 0.7894784210526316 0.7277034408713002 -1.6882833234244718 +1.00000001 0.7894784210526316 0.727703453398357 -1.688283343779765 +1.99999999 0.7894784210526316 4.071058568734327 -7.977813521160524 +2.0 0.7894784210526316 4.071058642376354 -7.977813667062976 +2.00000001 0.7894784210526316 4.071058716018381 -7.977813812965431 +2.99999999 0.7894784210526316 28.452786844330806 -60.091747774245306 +3.0 0.7894784210526316 28.452787436240953 -60.09174908712093 +3.00000001 0.7894784210526316 28.452788028151115 -60.091750399996585 +-1.0e-8 0.8684252631578947 -4.158803707010707e-9 -0.758601941360789 +1.0e-8 0.8684252631578947 4.158803707010707e-9 -0.758601941360789 +0.99999999 0.8684252631578947 0.5872819084776745 -1.383547833089991 +1.0 0.8684252631578947 0.5872819183316489 -1.3835478485879724 +1.00000001 0.8684252631578947 0.5872819281856234 -1.383547864085954 +1.99999999 0.8684252631578947 3.099596230982321 -5.956712345615788 +2.0 0.8684252631578947 3.0995962848854455 -5.9567124488634535 +2.00000001 0.8684252631578947 3.099596338788571 -5.956712552111121 +2.99999999 0.8684252631578947 20.06225114629942 -40.791740495677466 +3.0 0.8684252631578947 20.06225154720184 -40.791741348178505 +3.00000001 0.8684252631578947 20.062251948104265 -40.791742200679565 +-1.0e-8 0.9473721052631578 -3.463246204633292e-9 -0.6589920474143592 +1.0e-8 0.9473721052631578 3.463246204633292e-9 -0.6589920474143592 +0.99999999 0.9473721052631578 0.4794044550011185 -1.1497744148416618 +1.0 0.9473721052631578 0.4794044628618897 -1.1497744268580226 +1.00000001 0.9473721052631578 0.4794044707226609 -1.1497744388743836 +1.99999999 0.9473721052631578 2.403272135336257 -4.554767303366098 +2.0 0.9473721052631578 2.4032721756358426 -4.554767378374566 +2.00000001 0.9473721052631578 2.403272215935429 -4.554767453383036 +2.99999999 0.9473721052631578 14.522289239372824 -28.605321711788374 +3.0 0.9473721052631578 14.522289518857836 -28.605322284891624 +3.00000001 0.9473721052631578 14.522289798342852 -28.60532285799489 +-1.0e-8 1.026318947368421 -2.9051240120769176e-9 -0.5757367066850855 +1.0e-8 1.026318947368421 2.9051240120769176e-9 -0.5757367066850855 +0.99999999 1.026318947368421 0.395131080099377 -0.9665030122099901 +1.0 1.026318947368421 0.39513108644467426 -0.9665030216696985 +1.00000001 1.026318947368421 0.3951310927899715 -0.9665030311294068 +1.99999999 1.026318947368421 1.8919416422749398 -3.552346459353305 +2.0 1.026318947368421 1.8919416729451535 -3.5523465150561577 +2.00000001 1.026318947368421 1.8919417036153676 -3.552346570759011 +2.99999999 1.026318947368421 10.745439398886303 -20.61205796498971 +3.0 1.026318947368421 10.74543959850341 -20.612058361620896 +3.00000001 1.026318947368421 10.74543979812052 -20.612058758252086 +-1.0e-8 1.1052657894736841 -2.452344035634794e-9 -0.5054178274945208 +1.0e-8 1.1052657894736841 2.452344035634794e-9 -0.5054178274945208 +0.99999999 1.1052657894736841 0.3283640489068752 -0.8202112554633069 +1.0 1.1052657894736841 0.32836405408083497 -0.8202112630070297 +1.00000001 1.1052657894736841 0.32836405925479467 -0.8202112705507526 +1.99999999 1.1052657894736841 1.508744550712833 -2.8170662281283567 +2.0 1.1052657894736841 1.508744574411181 -2.817066270271266 +2.00000001 1.1052657894736841 1.5087445981095293 -2.817066312414176 +2.99999999 1.1052657894736841 8.100217810600931 -15.198256622790245 +3.0 1.1052657894736841 8.10021795614219 -15.198256904142767 +3.00000001 1.1052657894736841 8.100218101683451 -15.198257185495295 +-1.0e-8 1.1842126315789474 -2.0815492109678238e-9 -0.44550294924928646 +1.0e-8 1.1842126315789474 2.0815492109678238e-9 -0.44550294924928646 +0.99999999 1.1842126315789474 0.27482815614891276 -0.7016567946560259 +1.0 1.1842126315789474 0.27482816040471997 -0.7016568007388112 +1.00000001 1.1842126315789474 0.2748281646605272 -0.7016568068215964 +1.99999999 1.1842126315789474 1.2165570260953862 -2.2658831318073247 +2.0 1.1842126315789474 1.216557044647576 -2.265883164203767 +2.00000001 1.1842126315789474 1.216557063199766 -2.2658831966002104 +2.99999999 1.1842126315789474 6.204460312818431 -11.42970724247573 +3.0 1.1842126315789474 6.20446042083182 -11.429707446326844 +3.00000001 1.1842126315789474 6.204460528845212 -11.429707650177964 +-1.0e-8 1.2631594736842104 -1.7754031889777355e-9 -0.39407335384831066 +1.0e-8 1.2631594736842104 1.7754031889777355e-9 -0.39407335384831066 +0.99999999 1.2631594736842104 0.23145496436892718 -0.6043388789726469 +1.0 1.2631594736842104 0.2314549678962553 -0.60433888392474 +1.00000001 1.2631594736842104 0.23145497142358343 -0.6043388888768331 +1.99999999 1.2631594736842104 0.9904185900492432 -1.8449007997287565 +2.0 1.2631594736842104 0.9904186047389736 -1.8449008249790515 +2.00000001 1.2631594736842104 0.9904186194287041 -1.844900850229347 +2.99999999 1.2631594736842104 4.818604018890371 -8.74365800106262 +3.0 1.2631594736842104 4.8186041002983355 -8.743658151499659 +3.00000001 1.2631594736842104 4.818604181706303 -8.7436583019367 +-1.0e-8 1.3421063157894737 -1.520818566887326e-9 -0.34964720694443385 +1.0e-8 1.3421063157894737 1.520818566887326e-9 -0.34964720694443385 +0.99999999 1.3421063157894737 0.19599733177967094 -0.5235705007250242 +1.0 1.3421063157894737 0.19599733472286174 -0.5235705047906077 +1.00000001 1.3421063157894737 0.19599733766605257 -0.5235705088561914 +1.99999999 1.3421063157894737 0.8131167067232538 -1.5180967197442556 +2.0 1.3421063157894737 0.8131167184714848 -1.5180967396635003 +2.00000001 1.3421063157894737 0.813116730219716 -1.5180967595827453 +2.99999999 1.3421063157894737 3.787851594887129 -6.7893080054234085 +3.0 1.3421063157894737 3.7878516570786753 -6.789308118241521 +3.00000001 1.3421063157894737 3.787851719270223 -6.789308231059635 +-1.0e-8 1.4210531578947367 -1.3077687569951307e-9 -0.31106097679064526 +1.0e-8 1.4210531578947367 1.3077687569951307e-9 -0.31106097679064526 +0.99999999 1.4210531578947367 0.16678050999491892 -0.4558986509247457 +1.0 1.4210531578947367 0.16678051246533562 -0.4558986542873327 +1.00000001 1.4210531578947367 0.1667805149357523 -0.4558986576499197 +1.99999999 1.4210531578947367 0.6725173929857242 -1.2607697468355756 +2.0 1.4210531578947367 0.6725174024649351 -1.2607697627169807 +2.00000001 1.4210531578947367 0.6725174119441462 -1.260769778598386 +2.99999999 1.4210531578947367 3.0095004968727292 -5.341357345402829 +3.0 1.4210531578947367 3.0095005449554715 -5.341357431217093 +3.00000001 1.4210531578947367 3.0095005930382146 -5.341357517031358 +-1.0e-8 1.5 -1.1284720614977178e-9 -0.2773878004568438 +1.0e-8 1.5 1.1284720614977178e-9 -0.2773878004568438 +0.99999999 1.5 0.14253703968040865 -0.39873076015258935 +1.0 1.5 0.14253704176501716 -0.39873076295208826 +1.00000001 1.5 0.14253704384962568 -0.3987307657515872 +1.99999999 1.5 0.5598997819202102 -1.0555957386839712 +2.0 1.5 0.559899789629148 -1.0555957514657115 +2.00000001 1.5 0.5598997973380858 -1.0555957642474523 +2.99999999 1.5 2.4138110608791763 -4.251263302105232 +3.0 1.5 2.413811098451613 -4.25126336820581 +3.00000001 1.5 2.41381113602405 -4.251263434306388 From e317d3a3afdd2f00bccfcccdc8cb8d047c49fee5 Mon Sep 17 00:00:00 2001 From: Chris Geoga Date: Thu, 18 May 2023 13:48:11 -0400 Subject: [PATCH 15/17] Rolling back to the old power series. --- src/BesselFunctions/besselk.jl | 48 ++++++++++++++++++++++------------ 1 file changed, 31 insertions(+), 17 deletions(-) diff --git a/src/BesselFunctions/besselk.jl b/src/BesselFunctions/besselk.jl index 52ab073..58e81fc 100644 --- a/src/BesselFunctions/besselk.jl +++ b/src/BesselFunctions/besselk.jl @@ -500,24 +500,38 @@ besselk_power_series(v, x::ComplexF32) = ComplexF32(besselk_power_series(v, Comp function besselk_power_series(v, x::ComplexOrReal{T}) where T Math.isnearint(v) && return besselk_power_series_int(v, x) - MaxIter = 5000 - gam = gamma(v) - ngam = π / (sinpi(-abs(v)) * gam * v) - - s1, s2 = zero(T), zero(T) - t1, t2 = one(T), one(T) - - for k in 1:MaxIter - s1 += t1 - s2 += t2 - t1 *= x^2 / (4k * (k - v)) - t2 *= x^2 / (4k * (k + v)) - abs(t1) < eps(T) && break + MaxIter = 1000 + S = eltype(x) + v, x = S(v), S(x) + + z = x / 2 + zz = z * z + logz = log(z) + xd2_v = exp(v*logz) + xd2_nv = inv(xd2_v) + + # use the reflection identify to calculate gamma(-v) + # use relation gamma(v)*v = gamma(v+1) to avoid two gamma calls + gam_v = gamma(v) + gam_nv = π / (sinpi(-abs(v)) * gam_v * v) + gam_1mv = -gam_nv * v + gam_1mnv = gam_v * v + + _t1 = gam_v * xd2_nv * gam_1mv + _t2 = gam_nv * xd2_v * gam_1mnv + (xd2_pow, fact_k, out) = (one(S), one(S), zero(S)) + for k in 0:MaxIter + t1 = xd2_pow * T(0.5) + tmp = muladd(_t1, gam_1mnv, _t2 * gam_1mv) + tmp *= inv(gam_1mv * gam_1mnv * fact_k) + term = t1 * tmp + out += term + abs(term / out) < eps(T) && break + (gam_1mnv, gam_1mv) = (gam_1mnv*(one(S) + v + k), gam_1mv*(one(S) - v + k)) + xd2_pow *= zz + fact_k *= k + one(S) end - - xpv = (x/2)^v - s = gam * s1 + xpv^2 * ngam * s2 - return s / (2*xpv) + return out end besselk_power_series_cutoff(nu, x::Float64) = x < 2.0 || nu > 1.6x - 1.0 besselk_power_series_cutoff(nu, x::Float32) = x < 10.0f0 || nu > 1.65f0*x - 8.0f0 From 1578b6f0336a67c8c821f78989664cc035d714b1 Mon Sep 17 00:00:00 2001 From: Chris Geoga Date: Thu, 18 May 2023 13:57:44 -0400 Subject: [PATCH 16/17] Update README for test data generation --- test/data/besselk/enzyme/README.md | 42 ++++++++++++++++++++++-------- 1 file changed, 31 insertions(+), 11 deletions(-) diff --git a/test/data/besselk/enzyme/README.md b/test/data/besselk/enzyme/README.md index 4732d58..f038458 100644 --- a/test/data/besselk/enzyme/README.md +++ b/test/data/besselk/enzyme/README.md @@ -1,20 +1,19 @@ ```julia -using FiniteDifferences, ArbNumerics, DelimitedFiles, SpecialFunctions - -if !(@isdefined vx) - vgrid = sort(vcat(range(0.05, 15.0, length=20), [0.5, 1.5, 2.5, 3.5])) - xgrid = range(15.0, 30.0, length=20) - const vx = vec(collect(Iterators.product(vgrid, xgrid))) -end +using ArbNumerics, DelimitedFiles, SpecialFunctions +# Because you can't just use FiniteDifferences due to the "numerical noise". simplefd(f,x,h=ArbReal(1e-40)) = (f(x+h)-f(x))/h ArbNumerics.setprecision(ArbReal, digits=100) +arb_besselk(v,x) = ArbNumerics.besselk(ArbReal(v), ArbReal(x)) +arb_besselkx(v,x) = arb_besselk(v,x)*ArbNumerics.exp(ArbReal(x)) -function arb_besselkx(v,x) - (av, ax) = ArbReal.((v,x)) - ArbNumerics.besselk(av, ax)*ArbNumerics.exp(ax) -end +# +# besselkx_levin test: +# +vgrid = sort(vcat(range(0.05, 15.0, length=20), [0.5, 1.5, 2.5, 3.5])) +xgrid = range(15.0, 30.0, length=20) +vx = vec(collect(Iterators.product(vgrid, xgrid))) ref_values = map(vx) do vxj (v,x) = vxj @@ -29,4 +28,25 @@ out_matrix = hcat(getindex.(vx, 1), # test v argument getindex.(ref_values, 2)) # test d/dx value writedlm("besselkx_levin_enzyme_tests.csv", out_matrix) + +# +# besselk_power_series test: +# +vgrid = [-1e-8, 1e-8, 1-1e-8, 1.0, 1+1e-8, 2-1e-8, 2.0, 2+1e-8, 3-1e-8, 3.0, 3+1e-8] +xgrid = range(1e-5, 1.5, length=20) +vx = vec(collect(Iterators.product(vgrid, xgrid))) + +ref_values = map(vx) do vxj + (v,x) = vxj + dx = simplefd(_x->arb_besselk(v, _x), x) # NOT besselkx! + dv = simplefd(_v->arb_besselk(_v, x), v) # NOT besselkx! + Float64.((dv, dx)) +end + +out_matrix = hcat(getindex.(vx, 1), # test v argument + getindex.(vx, 2), # test x argument + getindex.(ref_values, 1), # test d/dv value + getindex.(ref_values, 2)) # test d/dx value + +writedlm("besselk_power_series_enzyme_tests.csv", out_matrix) ``` From 27d257fb8a3f51c411bb1bbb7cc1a43ee66c2233 Mon Sep 17 00:00:00 2001 From: Chris Geoga Date: Thu, 18 May 2023 14:11:49 -0400 Subject: [PATCH 17/17] Fix manual `sinpi` rule, but CI and compat to 1.9. --- .github/workflows/CI.yml | 2 +- Project.toml | 2 +- ext/BesselsEnzymeCoreExt.jl | 21 ++------------------- src/BesselFunctions/besselk.jl | 3 ++- 4 files changed, 6 insertions(+), 22 deletions(-) diff --git a/.github/workflows/CI.yml b/.github/workflows/CI.yml index 9463bdc..8223662 100644 --- a/.github/workflows/CI.yml +++ b/.github/workflows/CI.yml @@ -21,7 +21,7 @@ jobs: matrix: version: - '1' - - '1.8' + - '1.9' - 'nightly' os: - ubuntu-latest diff --git a/Project.toml b/Project.toml index edd7ab9..99624ca 100644 --- a/Project.toml +++ b/Project.toml @@ -7,7 +7,7 @@ SIMDMath = "5443be0b-e40a-4f70-a07e-dcd652efc383" [compat] SIMDMath = "0.2.5" -julia = "1.8" +julia = "1.9" [extras] Test = "8dfed614-e22c-5e08-85e1-65c5234f0b40" diff --git a/ext/BesselsEnzymeCoreExt.jl b/ext/BesselsEnzymeCoreExt.jl index 74f87b1..7aae17f 100644 --- a/ext/BesselsEnzymeCoreExt.jl +++ b/ext/BesselsEnzymeCoreExt.jl @@ -1,7 +1,5 @@ module BesselsEnzymeCoreExt - # TODO (cg 2023/05/08 10:02): Compat of any kind. - using Bessels, EnzymeCore using EnzymeCore.EnzymeRules using Bessels.Math @@ -9,22 +7,6 @@ module BesselsEnzymeCoreExt # A manual method that separately transforms the `val` and `dval`, because # sometimes the `val` can converge while the `dval` hasn't, so just using an # early return or something can give incorrect derivatives in edge cases. - # - # https://github.com/JuliaMath/Bessels.jl/issues/96 - # - # and links with for discussion. - # - # TODO (cg 2023/05/08 10:00): I'm not entirely sure how best to "generalize" - # this to cases like a return type of DuplicatedNoNeed, or something being a - # `Enzyme.Const`. These shouldn't in principle affect the "point" of this - # function (which is just to check for convergence before applying a - # function), but on its face this approach would mean I need a lot of - # hand-written extra methods. I have an open issue on the Enzyme.jl repo at - # - # https://github.com/EnzymeAD/Enzyme.jl/issues/786 - # - # that gets at this problem a bit. But it's a weird request and I'm sure Billy - # has a lot of asks on his time. function EnzymeRules.forward(func::Const{typeof(levin_transform)}, ::Type{<:Duplicated}, s::Duplicated, @@ -39,7 +21,8 @@ module BesselsEnzymeCoreExt function EnzymeRules.forward(func::Const{typeof(sinpi)}, ::Type{<:Duplicated}, x::Duplicated) - Duplicated(sinpi(x.val), pi*cospi(x.val)) + (sp, cp) = sincospi(x.val) + Duplicated(sp, pi*cp*x.dval) end function EnzymeRules.forward(func::Const{typeof(sinpi)}, diff --git a/src/BesselFunctions/besselk.jl b/src/BesselFunctions/besselk.jl index 58e81fc..4bed7e4 100644 --- a/src/BesselFunctions/besselk.jl +++ b/src/BesselFunctions/besselk.jl @@ -513,7 +513,8 @@ function besselk_power_series(v, x::ComplexOrReal{T}) where T # use the reflection identify to calculate gamma(-v) # use relation gamma(v)*v = gamma(v+1) to avoid two gamma calls gam_v = gamma(v) - gam_nv = π / (sinpi(-abs(v)) * gam_v * v) + #gam_nv = π / (sin(-pi*abs(v)) * gam_v * v) # not using sinpi here to avoid Enzyme bug + gam_nv = π / (sinpi(-abs(v)) * gam_v * v) # not using sinpi here to avoid Enzyme bug gam_1mv = -gam_nv * v gam_1mnv = gam_v * v