merging into main

Author: shahriariravanian

Project.toml

@@ -1,14 +1,13 @@
 name = "SymbolicNumericIntegration"
 uuid = "78aadeae-fbc0-11eb-17b6-c7ec0477ba9e"
 authors = ["Shahriar Iravanian <[email protected]>"]
-version = "1.2.1"
+version = "1.3.0"
 
 [deps]
 DataDrivenDiffEq = "2445eb08-9709-466a-b3fc-47e12bd697a2"
 DataDrivenSparse = "5b588203-7d8b-4fab-a537-c31a7f73f46b"
 DataStructures = "864edb3b-99cc-5e75-8d2d-829cb0a9cfe8"
 LinearAlgebra = "37e2e46d-f89d-539d-b4ee-838fcccc9c8e"
-Optim = "429524aa-4258-5aef-a3af-852621145aeb"
 SpecialFunctions = "276daf66-3868-5448-9aa4-cd146d93841b"
 Statistics = "10745b16-79ce-11e8-11f9-7d13ad32a3b2"
 SymbolicUtils = "d1185830-fcd6-423d-90d6-eec64667417b"
@@ -18,7 +17,6 @@ Symbolics = "0c5d862f-8b57-4792-8d23-62f2024744c7"
 DataDrivenDiffEq = "1.0"
 DataDrivenSparse = "0.1.2"
 DataStructures = "0.18"
-Optim = "1"
 SymbolicUtils = "1"
 Symbolics = "5"
 SpecialFunctions = "2"

src/SymbolicNumericIntegration.jl

@@ -33,7 +33,7 @@ export Ei, Si, Ci, Li
 
 include("numeric_utils.jl")
 include("sparse.jl")
-include("optim.jl")
+# include("optim.jl")
 include("integral.jl")
 
 export integrate, generate_basis

src/integral.jl

@@ -43,8 +43,8 @@ Keyword Arguments:
 - `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*)
+- `homotopy` (default: `true`): *deprecated*, will be removed in a future version
+- `use_optim` (default: `false`): *deprecated*, will be removed in a future version
 - `detailed` (default: `true`): `(solved, unsolved, err)` output format. If `detailed=false`, only the final integral is returned. 
 
 Output:
@@ -67,7 +67,9 @@ function integrate(eq, x = nothing;
     complex_plane = true,
     homotopy = true, 
     use_optim = false, 
-    detailed = true)
+    detailed = true)    
+    
+    deprecation_warnings(; homotopy, use_optim)   
 
     eq = expand(eq)
 
@@ -309,16 +311,13 @@ end
 
 #################################################################################
 
-"""
-	integrate_basis(eq, x; kwargs...)
-	
-is used for debugging and should not be called in the course of normal execution
-"""
-function integrate_basis(eq, x = var(eq); plan = default_plan())
-    eq = cache(expand(eq))
-    basis = generate_basis(eq, x, false)
-    n = length(basis)
-    A, X = init_basis_matrix(eq, x, basis; plan)
-    return basis, A, X
+function deprecation_warnings(; use_optim=false, homotopy=true)
+    if use_optim
+        @warn("use_optim is deprecated and will be removed in a future version")
+    end
+    
+    if !homotopy
+        @warn("homotopy is deprecated and will be removed in a future version")
+    end
 end
 

src/optim.jl

@@ -199,7 +199,8 @@ end
 
 function expand_basis_symbolic(basis, x)
     b = sum(basis)
-    basis = split_terms(expand((1+x)*(b + diff(b, x))), x)    
+    Δb = sum(diff(y, x) for y in basis)
+    basis = split_terms(expand((1+x)*(b + Δb)), x)    
 
     return basis
 end
@@ -256,7 +257,11 @@ function Problem(eq, x; plan = default_plan(), num_steps=1)
     ker = best_hints(eq, x, basis; plan)
 
     if ker == nothing
-        return nothing
+        basis = expand_basis_symbolic(basis, x)
+        ker = best_hints(eq, x, basis; plan)
+        if ker == nothing
+            return nothing
+        end
     end
 
     ker = [atomize(y, x)[2] for y in ker]

src/symbolic.jl

@@ -197,6 +197,11 @@ basic_integrals = [
     1 / log(3x - 1),
     1 / (x * log(log(x))),
     x / log(x^2),
+        # from Geddes & Stefanus
+    exp(1 / (1 + log(x)))*((2*log(x)+1)/x),
+    1 / (1 + exp(x)),
+    sin(2x) * exp(x),
+    -10*exp(x^-10) / x^11,
     # bypass = true
     β,      # turn of bypass = true
     (log(x - 1) + (x - 1)^-1) * log(x),
@@ -229,7 +234,7 @@ sym_integrals = [
     log(a * x)^2,
     log(a + x),
     log(a + x^2) * x,    
-    x^2 * log(a * x + b)^2,
+    # x^2 * log(a * x + b)^2,
     exp(a * x),
     x * exp(a * x),
     x^2 * exp(a * x),
@@ -258,7 +263,7 @@ sym_integrals = [
     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),