From 15fe5759799f47e4b3c04f07f18952af607e79f1 Mon Sep 17 00:00:00 2001 From: Hendrik Ranocha Date: Tue, 16 Aug 2022 19:24:16 +0200 Subject: [PATCH] symbolic division by h and latexify (#72) * test symbolic division by h and latexify * symbolic division by h in tutorial * set version to v0.1.30 * format --- Project.toml | 2 +- docs/src/tutorials/symbolic_computations.md | 42 +++++++++- test/runtests.jl | 87 ++++++++++++++++----- 3 files changed, 108 insertions(+), 23 deletions(-) diff --git a/Project.toml b/Project.toml index 4c56c5fa..d935f9b0 100644 --- a/Project.toml +++ b/Project.toml @@ -1,7 +1,7 @@ name = "BSeries" uuid = "ebb8d67c-85b4-416c-b05f-5f409e808f32" authors = ["Hendrik Ranocha and contributors"] -version = "0.1.29" +version = "0.1.30" [deps] Latexify = "23fbe1c1-3f47-55db-b15f-69d7ec21a316" diff --git a/docs/src/tutorials/symbolic_computations.md b/docs/src/tutorials/symbolic_computations.md index aebae6f5..1bca4979 100644 --- a/docs/src/tutorials/symbolic_computations.md +++ b/docs/src/tutorials/symbolic_computations.md @@ -69,6 +69,22 @@ using Latexify latexify(series, reduce_order_by=1, dt=SymEngine.symbols("h"), cdot=false) |> println ``` +The argument `dt=SymEngine.symbols("h")` is not required here. + +```@example modified-equation-symengine +using Latexify +latexify(series, reduce_order_by=1, cdot=false) |> println +``` + +We can obtain basically the same result (up to other notations of the coefficients) +by dividing the B-series `series` by a symbolic variable for the time step size +as follows. + +```@example modified-equation-symengine +using Latexify +h = SymEngine.symbols("h") +latexify(series / h, cdot=false) |> println +``` We can also use other packages for the symbolic computations, of course. SymPy.jl often provides very clean expressions. @@ -92,13 +108,25 @@ using Latexify latexify(series, reduce_order_by=1, dt=SymPy.symbols("h"), cdot=false) |> println ``` +We can also use the simplified versions. + +```@example modified-equation-sympy +using Latexify +latexify(series, reduce_order_by=1, cdot=false) |> println +``` + +```@example modified-equation-sympy +using Latexify +h = SymPy.symbols("h", real=true) +latexify(series / h, cdot=false) |> println +``` Alternatively, we can also use Symbolics.jl. ```@example modified-equation-symbolics using BSeries, Symbolics -Symbolics.@variables α h +Symbolics.@variables α A = [0 0; 1/(2α) 0] b = [1-α, α] c = [0, 1/(2α)] @@ -106,14 +134,20 @@ c = [0, 1/(2α)] series = modified_equation(A, b, c, 3) ``` +```@example modified-equation-symbolics +using Latexify +Symbolics.@variables h +latexify(series / h, cdot=false) |> println +``` + ## Working with rational coefficients High-order Runge-Kutta methods are often given in terms of rational coefficients that have large denominators. -It's often best to use the rational coefficients in calculations, rather than converting to floating-point, since otherwise +It's often best to use the rational coefficients in calculations, rather than converting to floating-point, since otherwise terms that ought to cancel may not cancel exactly. Rational coefficients can be entered by using `//` for the division operation. -By default, Julia uses 64-bit integers for the numerator and denominator of a rational on 64-bit systems. In practical -calculations with high-order RK methods, the denominators may become too large to be represented with 64 bits, +By default, Julia uses 64-bit integers for the numerator and denominator of a rational on 64-bit systems. In practical +calculations with high-order RK methods, the denominators may become too large to be represented with 64 bits, leading to overflow. This can be remedied by specifying higher precision when entering the coefficients: ```@example int128-coefficients-symbolics diff --git a/test/runtests.jl b/test/runtests.jl index 474c8cde..0d034f25 100644 --- a/test/runtests.jl +++ b/test/runtests.jl @@ -31,30 +31,81 @@ using Aqua: Aqua @test_nowarn latexify(series_integrator, dt = SymEngine.symbols("h")) @testset "SymEngine.jl" begin - α = SymEngine.symbols("α") - A = [0 0; 1/(2 * α) 0] - b = [1 - α, α] - c = [0, 1 / (2 * α)] - series_integrator = bseries(A, b, c, 3) - @test_nowarn latexify(series_integrator) + @testset "Symbolic coefficients" begin + α = SymEngine.symbols("α") + A = [0 0; 1/(2 * α) 0] + b = [1 - α, α] + c = [0, 1 / (2 * α)] + series_integrator = bseries(A, b, c, 3) + @test_nowarn latexify(series_integrator) + end + + @testset "Divide by h" begin + A = [0 0; 1//2 0] + b = [0, 1] + c = [0, 1 // 2] + series_integrator = bseries(A, b, c, 1) + @test_nowarn latexify(series_integrator) + + h = SymEngine.symbols("h") + coefficients = modified_equation(series_integrator) + val1 = @test_nowarn latexify(coefficients, reduce_order_by = 1, + cdot = false) + val2 = @test_nowarn latexify(coefficients / h, cdot = false) + @test val1 == val2 + end end @testset "SymPy.jl" begin - α = SymPy.symbols("α", real = true) - A = [0 0; 1/(2 * α) 0] - b = [1 - α, α] - c = [0, 1 / (2 * α)] - series_integrator = bseries(A, b, c, 3) - @test_nowarn latexify(series_integrator) + @testset "Symbolic coefficients" begin + α = SymPy.symbols("α", real = true) + A = [0 0; 1/(2 * α) 0] + b = [1 - α, α] + c = [0, 1 / (2 * α)] + series_integrator = bseries(A, b, c, 3) + @test_nowarn latexify(series_integrator) + end + + @testset "Divide by h" begin + A = [0 0; 1//2 0] + b = [0, 1] + c = [0, 1 // 2] + series_integrator = bseries(A, b, c, 1) + @test_nowarn latexify(series_integrator) + + h = SymPy.symbols("h", real = true) + coefficients = modified_equation(series_integrator) + val1 = @test_nowarn latexify(coefficients, reduce_order_by = 1, + cdot = false) + val2 = @test_nowarn latexify(coefficients / h, cdot = false) + @test val1 == val2 + end end @testset "Symbolics.jl" begin - Symbolics.@variables α - A = [0 0; 1/(2 * α) 0] - b = [1 - α, α] - c = [0, 1 / (2 * α)] - series_integrator = bseries(A, b, c, 3) - @test_nowarn latexify(series_integrator) + @testset "Symbolic coefficients" begin + Symbolics.@variables α + A = [0 0; 1/(2 * α) 0] + b = [1 - α, α] + c = [0, 1 / (2 * α)] + series_integrator = bseries(A, b, c, 3) + @test_nowarn latexify(series_integrator) + end + + @testset "Divide by h" begin + A = [0 0; 1//2 0] + b = [0, 1] + c = [0, 1 // 2] + series_integrator = bseries(A, b, c, 1) + @test_nowarn latexify(series_integrator) + + Symbolics.@variables h + coefficients = modified_equation(series_integrator) + val1 = @test_nowarn latexify(coefficients, reduce_order_by = 1, + cdot = false) + val2 = @test_nowarn latexify(coefficients / h, cdot = false) + @test val1 == val2 + end end end # @testset "latexify"