merged with hyint for symbolic integration

Author: shahriariravanian

Project.toml

@@ -1,7 +1,7 @@
 name = "SymbolicNumericIntegration"
 uuid = "78aadeae-fbc0-11eb-17b6-c7ec0477ba9e"
 authors = ["Shahriar Iravanian <[email protected]>"]
-version = "1.1.0"
+version = "1.2.0"
 
 [deps]
 DataDrivenDiffEq = "2445eb08-9709-466a-b3fc-47e12bd697a2"

README.md

@@ -24,13 +24,16 @@ the documentation, which contains the unreleased features.
 using Symbolics
 using SymbolicNumericIntegration
 
-@variables x
+@variables x a b
 
-integrate(3x^3 + 2x - 5)
-```
+julia> integrate(3x^3 + 2x - 5)
+(x^2 + (3//4)*(x^4) - (5//1)*x, 0, 0)
 
-```
-(x^2 + (3//4)*(x^4) - (5x), 0, 0)
+julia> integrate(exp(a*x), x; symbolic=true)
+(exp(a*x) / a, 0, 0)
+
+julia> integrate(sin(a*x) * cos(b*x), x; symbolic=true, detailed=false)
+(-a*cos(a*x)*cos(b*x) - b*sin(a*x)*sin(b*x)) / (a^2 - (b^2))
 ```
 
 # Citation

docs/make.jl

@@ -9,14 +9,7 @@ makedocs(sitename = "SymbolicNumericIntegration.jl",
          authors = "Shahriar Iravanian",
          modules = [SymbolicNumericIntegration],
          clean = true, doctest = false, linkcheck = true,
-         strict = [
-             :doctest,
-             :linkcheck,
-             :parse_error,
-             :example_block,
-             # Other available options are
-             # :autodocs_block, :cross_references, :docs_block, :eval_block, :example_block, :footnote, :meta_block, :missing_docs, :setup_block
-         ],
+         warnonly = true,
          format = Documenter.HTML(analytics = "UA-90474609-3",
                                   assets = ["assets/favicon.ico"],
                                   canonical = "https://docs.sciml.ai/SymbolicNumericIntegration/stable/"),

docs/src/index.md

@@ -2,10 +2,13 @@
 
 **SymbolicNumericIntegration.jl** is a hybrid symbolic/numerical integration package that works on the [Julia Symbolics](https://docs.sciml.ai/Symbolics/stable/) expressions.
 
-**SymbolicNumericIntegration.jl** uses a randomized algorithm based on a hybrid of the *method of undetermined coefficients* and *sparse regression* and can solve a large subset of basic standard integrals (polynomials, exponential/logarithmic, trigonometric and hyperbolic, inverse trigonometric and hyperbolic, rational and square root).
-The basis of how it works and the theory of integration using the Symbolic-Numeric methods refer to [Basis of Symbolic-Numeric Integration](https://github.com/SciML/SymbolicNumericIntegration.jl/blob/main/docs/theory.ipynb).
+**SymbolicNumericIntegration.jl** uses a randomized algorithm based on a hybrid of the *method of undetermined coefficients* and *sparse regression* and can solve a large subset of basic standard integrals (polynomials, exponential/logarithmic, trigonometric and hyperbolic, inverse trigonometric and hyperbolic, rational and square root). The symbolic part of the algorithm is similar (but not identical) to Risch-Bronstein's poor man's integrator and generates a list of ansatzes (candidate terms). The numerical part uses sparse regression adopted from the Sparse identification of nonlinear dynamics (SINDy) algorithm to prune down the ansatzes and find the corresponding coefficients. The basis of how it works and the theory of integration using the Symbolic-Numeric methods refer to [Basis of Symbolic-Numeric Integration](https://github.com/SciML/SymbolicNumericIntegration.jl/blob/main/docs/theory.ipynb). 
 
-Function `integrate` returns the integral of a univariate expression with *constant* real or complex coefficients. `integrate` returns a tuple with three values. The first one is the solved integral, the second one is the sum of the unsolved terms, and the third value is the residual error. If `integrate` is successful, the unsolved portion is reported as 0.
+[hyint](https://github.com/siravan/hyint) is the python counterpart of **SymbolicNumericIntegration.jl** that works with **sympy** expressions. 
+
+Originally, **SymbolicNumericIntegration.jl** expected only univariate expression with *constant* real or complex coefficients as input. As of version 1.2, `integrate` function accepts symbolic constants for a subset of solvable integrals. 
+
+`integrate` returns a tuple with three values. The first one is the solved integral, the second one is the sum of the unsolved terms, and the third value is the residual error. If `integrate` is successful, the unsolved portion is reported as 0. If we pass `detailed=false` to `integrate', the output is simplified to only the resulting integrals. In this case, `nothing` is returned if the integral cannot be solved. Note that the simplified output will become the default in a future version. 
 
 ## Installation
 
@@ -22,7 +25,7 @@ Examples:
 using Symbolics
 using SymbolicNumericIntegration
 
-@variables x
+@variables x a b
 ```
 
 ```julia
@@ -32,44 +35,31 @@ julia> integrate(3x^3 + 2x - 5)
 julia> integrate((5 + 2x)^-1)
 ((1//2)*log((5//2) + x), 0, 0.0)
 
-julia> integrate(1 / (6 + x^2 - (5x)))
-(log(x - 3) - log(x - 2), 0, 3.339372764128952e-16)
-
-julia> integrate(1 / (x^2 - 16))
-((1//8)*log(x - 4) - ((1//8)*log(4 + x)), 0, 1.546926788028958e-16)
-
-julia> integrate(x^2 / (16 + x^2))
-(x + 4atan((-1//4)*x), 0, 1.3318788420751984e-16)
-
-julia> integrate(x^2 / sqrt(4 + x^2))
-((1//2)*x*((4 + x^2)^0.5) - ((2//1)*log(x + sqrt(4 + x^2))), 0, 8.702422633074313e-17)
-
-julia> integrate(x^2 * log(x))
-((1//3)*log(x)*(x^3) - ((1//9)*(x^3)), 0, 0)
+# detailed simplifies the output to just the resulting integral
+julia> integrate(x^2 / (16 + x^2); detailed=false)
+x + 4atan((-1//4)*x)
 
-julia> integrate(x^2 * exp(x))
-(2exp(x) + exp(x)*(x^2) - (2x*exp(x)), 0, 0)
+julia> integrate(x^2 * log(x); detailed=false)
+(1//3)*(x^3)*log(x) - (1//9)*(x^3)
 
-julia> integrate(tan(2x))
-((-1//2)*log(cos(2x)), 0, 0)
+julia> integrate(sec(x) * tan(x); detailed=false)
+sec(x)
 
-julia> integrate(sec(x) * tan(x))
-(cos(x)^-1, 0, 0)
+# Here, a is a symbolic constant; therefore, we need to explicitly
+# define the independent variable (say, x). Also, we set
+# `symbolic=true` to force using the symbolic solver 
+julia> integrate(sin(a*x), x; detailed=false, symbolic=true)
+(-cos(a*x)) / a
 
-julia> integrate(cosh(2x) * exp(x))
-((2//3)*exp(x)*sinh(2x) - ((1//3)*exp(x)*cosh(2x)), 0, 7.073930088880992e-8)
+julia> integrate(x^2 * cos(a*x), x; detailed=false, symbolic=true)
+((a^2)*(x^2)*sin(a*x) + 2.0a*x*cos(a*x) - 2.0sin(a*x)) / (a^3)
 
-julia> integrate(cosh(x) * sin(x))
-((1//2)*sin(x)*sinh(x) - ((1//2)*cos(x)*cosh(x)), 0, 4.8956233716268386e-17)
+# multiple symbolic constants
+julia> integrate(cosh(a*x) * exp(b*x), x; detailed=false, symbolic=true)
+(a*sinh(a*x)*exp(b*x) - b*cosh(a*x)*exp(b*x)) / (a^2 - (b^2))
 
-julia> integrate(cosh(2x) * sin(3x))
-(0.153845sinh(2x)*sin(3x) - (0.23077cosh(2x)*cos(3x)), 0, 4.9807620877373405e-6)
-
-julia> integrate(log(log(x)) * (x^-1))
-(log(x)*log(log(x)) - log(x), 0, 0)
-
-julia> integrate(exp(x^2))
-(0, exp(x^2), Inf)    # as expected!
+julia> integrate(log(log(a*x)) / x, x; detailed=false, symbolic=true)
+log(a*x)*log(log(a*x)) - log(a*x)
 ```
 
 SymbolicNumericIntegration.jl exports some special integral functions (defined over Complex numbers) and uses them in solving integrals:
@@ -83,13 +73,13 @@ For examples:
 
 ```
 julia> integrate(exp(x + 1) / (x + 1))
-(SymbolicNumericIntegration.Ei(1 + x), 0, 1.1796119636642288e-16)
+(SymbolicNumericIntegration.Ei(1 + x), 0, 0.0)
 
-julia> integrate(x * cos(x^2 - 1) / (x^2 - 1))
-((1//2)*SymbolicNumericIntegration.Ci(x^2 - 1), 0, 2.7755575615628914e-17)
+julia> integrate(x * cos(a*x^2 - 1) / (a*x^2 - 1), x; detailed=false, symbolic=true)
+((1//2)*SymbolicNumericIntegration.Ci(a*(x^2) - 1)) / a
 
-julia> integrate(1 / (x*log(log(x))))
-(SymbolicNumericIntegration.Li(log(x)), 0, 1.1102230246251565e-16)
+julia> integrate(1 / (x*log(log(x))), x; detailed=false, symbolic=true)
+SymbolicNumericIntegration.Li(log(x))
 ```
 
 ```@docs
@@ -198,33 +188,3 @@ Pkg.status(; mode = PKGMODE_MANIFEST) # hide
 ```@raw html
 </details>
 ```
-
-```@raw html
-You can also download the 
-<a href="
-```
-
-```@eval
-using TOML
-version = TOML.parse(read("../../Project.toml", String))["version"]
-name = TOML.parse(read("../../Project.toml", String))["name"]
-link = "https://github.com/SciML/" * name * ".jl/tree/gh-pages/v" * version *
-       "/assets/Manifest.toml"
-```
-
-```@raw html
-">manifest</a> file and the
-<a href="
-```
-
-```@eval
-using TOML
-version = TOML.parse(read("../../Project.toml", String))["version"]
-name = TOML.parse(read("../../Project.toml", String))["name"]
-link = "https://github.com/SciML/" * name * ".jl/tree/gh-pages/v" * version *
-       "/assets/Project.toml"
-```
-
-```@raw html
-">project</a> file.
-```

src/integral.jl

@@ -29,14 +29,15 @@ Keyword Arguments:
 - `num_trials` (default: `10`): the number of trials in each step (no changes to the basis)
 - `show_basis` (default: `false`): if true, the basis (list of candidate terms) is printed
 - `bypass` (default: `false`): if true do not integrate terms separately but consider all at once
-- `symbolic` (default: `false`): try symbolic integration first (will be forced if `eq` has constant parameters)
+- `symbolic` (default: `false`): try symbolic integration first (will be forced if `eq` has symbolic constants)
 - `max_basis` (default: `100`): the maximum number of candidate terms to consider
 - `verbose` (default: `false`): print a detailed report
 - `complex_plane` (default: `true`): generate random test points on the complex plane (if false, the points will be on real axis)
 - `radius` (default: `1.0`): the radius of the disk in the complex plane to generate random test points
 - `opt` (default: `STLSQ(exp.(-10:1:0))`): the sparse regression optimizer (from DataDrivenSparse)
 - `homotopy` (default: `true`): use the homotopy algorithm to generate the basis (*deprecated*, will be removed in a future version)
 - `use_optim` (default: `false`): use Optim.jl `minimize` function instead of the STLSQ algorithm (*experimental*)
+- `detailed` (default: `true`): `(solved, unsolved, err)` output format. If `detailed=false`, only the final integral is returned. 
 
 Output:
 -------
@@ -52,11 +53,12 @@ function integrate(eq, x = nothing; abstol = 1e-6, num_steps = 2, num_trials = 1
     eq = expand(eq)
 
     if x == nothing
-        if !isempty(sym_consts(eq, x))
+        vars = get_variables(eq)
+        if length(vars) > 1
             error("Multiple symbolic variables detect. Please pass the independent variable to `integrate`")            
-        end
-        x = var(eq)
-        if x == nothing
+        elseif length(vars) == 1
+            x = vars[1]
+        else
             @syms 𝑥
             x = 𝑥
         end

src/sparse.jl

@@ -198,7 +198,7 @@ function hints(eq, x, basis; radius=5.0, abstol=1e-6, opt=STLSQ(exp.(-10:1:0)),
         end
         return h, err
     catch e
-        println("Error from hints", e)        
+        # println("Error from hints: ", e)        
     end
 
     return 0, Inf
@@ -207,7 +207,7 @@ end
 
 function best_hints(eq, x, basis; radius=5.0, abstol=1e-6, opt=STLSQ(exp.(-10:1:0)), complex_plane=true, num_trials=10)        
     H = []
-    L = Int[]
+    L = Int[]    
 
     for _ in 1:num_trials
         try 

test/runtests.jl

@@ -11,7 +11,7 @@ include("axiom.jl")
 
 ##############################################################################
 
-@variables x a b β
+@variables x a b c β
 
 """
     a list of basic standard integral tests
@@ -209,6 +209,7 @@ basic_integrals = [
 sym_integrals = [
     # Basic Forms
     a*x^2,
+    a*x + b*x^2 - c*x^3,
     a / x,
     1 / (a*x + 5),
     1 / (x + a)^2,
@@ -238,20 +239,23 @@ sym_integrals = [
     sec(a*x),
     x * cos(a*x),
     x^2 * cos(a*x),
-    sin(a*x)^2 * cos(b*x)^3,
     exp(a*x) * sin(b*x),
     x * exp(a*x) * sin(a*x),
     x * exp(a*x) * cos(b*x),
     cosh(a*x),
     exp(a*x) * cosh(b*x),
     cos(a*x) * cosh(b*x),
+    sin(a*x) * cos(b*x) * exp(c*x),
+    sin(a*x) * sinh(b*x) * exp(c*x),
     sec(a*x)^2 * tan(a*x),
     exp(a*x) / (1 + exp(a*x)),
+    exp(a*x) / exp(b*x),
     cos(exp(a*x)) * sin(exp(a*x)) * exp(a*x),
     1 / (x * log(a*x)),
     log(log(a*x)) / x,
     sin(log(a*x)),
-    x / (exp(a*x) - b),
+    x / (exp(a*x) - b),    
+    exp(a*x) / (b*exp(a*x) + c),
     exp(a*x) * exp(exp(a*x)),
     log(cos(a*x)) * tan(a*x),
     1 / (x^3 + a),
@@ -266,14 +270,14 @@ sym_integrals = [
 ]
 
 
-function test_integrals(basis=true, subs=nothing;  kw...)
+function test_integrals(basic=true, subs=nothing;  kw...)
     args = isempty(kw) ? Dict() : Dict(kw)
     args[:detailed] = false
     misses = []
     k = 1        
 
-    integrals = basis ? basic_integrals : sym_integrals    
-    args[:symbolic] = !basis
+    integrals = basic ? basic_integrals : sym_integrals    
+    args[:symbolic] = !basic
 
     for (i, eq) in enumerate(integrals)
         if isequal(eq, β)