Skip to content

Commit c9311c4

Browse files
formatted
1 parent e2b153a commit c9311c4

File tree

11 files changed

+379
-400
lines changed

11 files changed

+379
-400
lines changed

README.md

Lines changed: 3 additions & 5 deletions
Original file line numberDiff line numberDiff line change
@@ -21,19 +21,17 @@ the documentation, which contains the unreleased features.
2121
## Example
2222

2323
```julia
24+
julia> integrate(3x^3 + 2x - 5)
2425
using Symbolics
2526
using SymbolicNumericIntegration
2627

2728
@variables x a b
2829

29-
julia> integrate(3x^3 + 2x - 5)
30+
julia> integrate(exp(a * x), x; symbolic = true)
3031
(x^2 + (3//4)*(x^4) - (5//1)*x, 0, 0)
3132

32-
julia> integrate(exp(a*x), x; symbolic=true)
33+
julia> integrate(sin(a * x) * cos(b * x), x; symbolic = true, detailed = false)
3334
(exp(a*x) / a, 0, 0)
34-
35-
julia> integrate(sin(a*x) * cos(b*x), x; symbolic=true, detailed=false)
36-
(-a*cos(a*x)*cos(b*x) - b*sin(a*x)*sin(b*x)) / (a^2 - (b^2))
3735
```
3836

3937
# Citation

docs/make.jl

Lines changed: 9 additions & 9 deletions
Original file line numberDiff line numberDiff line change
@@ -6,14 +6,14 @@ using Documenter, SymbolicNumericIntegration
66
include("pages.jl")
77

88
makedocs(sitename = "SymbolicNumericIntegration.jl",
9-
authors = "Shahriar Iravanian",
10-
modules = [SymbolicNumericIntegration],
11-
clean = true, doctest = false, linkcheck = true,
12-
warnonly = true,
13-
format = Documenter.HTML(analytics = "UA-90474609-3",
14-
assets = ["assets/favicon.ico"],
15-
canonical = "https://docs.sciml.ai/SymbolicNumericIntegration/stable/"),
16-
pages = pages)
9+
authors = "Shahriar Iravanian",
10+
modules = [SymbolicNumericIntegration],
11+
clean = true, doctest = false, linkcheck = true,
12+
warnonly = true,
13+
format = Documenter.HTML(analytics = "UA-90474609-3",
14+
assets = ["assets/favicon.ico"],
15+
canonical = "https://docs.sciml.ai/SymbolicNumericIntegration/stable/"),
16+
pages = pages)
1717

1818
deploydocs(repo = "github.com/SciML/SymbolicNumericIntegration.jl.git";
19-
push_preview = true)
19+
push_preview = true)

docs/src/index.md

Lines changed: 15 additions & 12 deletions
Original file line numberDiff line numberDiff line change
@@ -2,13 +2,13 @@
22

33
**SymbolicNumericIntegration.jl** is a hybrid symbolic/numerical integration package that works on the [Julia Symbolics](https://docs.sciml.ai/Symbolics/stable/) expressions.
44

5-
**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).
5+
**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).
66

7-
[hyint](https://github.com/siravan/hyint) is the python counterpart of **SymbolicNumericIntegration.jl** that works with **sympy** expressions.
7+
[hyint](https://github.com/siravan/hyint) is the python counterpart of **SymbolicNumericIntegration.jl** that works with **sympy** expressions.
88

9-
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.
9+
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.
1010

11-
`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.
11+
`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.
1212

1313
## Installation
1414

@@ -36,29 +36,32 @@ julia> integrate((5 + 2x)^-1)
3636
((1//2)*log((5//2) + x), 0, 0.0)
3737

3838
# detailed simplifies the output to just the resulting integral
39-
julia> integrate(x^2 / (16 + x^2); detailed=false)
39+
40+
julia> integrate(x^2 / (16 + x^2); detailed = false)
4041
x + 4atan((-1//4)*x)
4142

42-
julia> integrate(x^2 * log(x); detailed=false)
43+
julia> integrate(x^2 * log(x); detailed = false)
4344
(1//3)*(x^3)*log(x) - (1//9)*(x^3)
4445

45-
julia> integrate(sec(x) * tan(x); detailed=false)
46+
julia> integrate(sec(x) * tan(x); detailed = false)
4647
sec(x)
4748

4849
# Here, a is a symbolic constant; therefore, we need to explicitly
4950
# define the independent variable (say, x). Also, we set
50-
# `symbolic=true` to force using the symbolic solver
51-
julia> integrate(sin(a*x), x; detailed=false, symbolic=true)
51+
# `symbolic=true` to force using the symbolic solver
52+
53+
julia> integrate(sin(a * x), x; detailed = false, symbolic = true)
5254
(-cos(a*x)) / a
5355

54-
julia> integrate(x^2 * cos(a*x), x; detailed=false, symbolic=true)
56+
julia> integrate(x^2 * cos(a * x), x; detailed = false, symbolic = true)
5557
((a^2)*(x^2)*sin(a*x) + 2.0a*x*cos(a*x) - 2.0sin(a*x)) / (a^3)
5658

5759
# multiple symbolic constants
58-
julia> integrate(cosh(a*x) * exp(b*x), x; detailed=false, symbolic=true)
60+
61+
julia> integrate(cosh(a * x) * exp(b * x), x; detailed = false, symbolic = true)
5962
(a*sinh(a*x)*exp(b*x) - b*cosh(a*x)*exp(b*x)) / (a^2 - (b^2))
6063

61-
julia> integrate(log(log(a*x)) / x, x; detailed=false, symbolic=true)
64+
julia> integrate(log(log(a * x)) / x, x; detailed = false, symbolic = true)
6265
log(a*x)*log(log(a*x)) - log(a*x)
6366
```
6467

src/homotopy.jl

Lines changed: 73 additions & 70 deletions
Original file line numberDiff line numberDiff line change
@@ -66,12 +66,13 @@ Symbolics.@register_symbolic Ci(z)
6666
Symbolics.@register_symbolic Li(z)
6767
Symbolics.@register_symbolic Erfi(z)
6868

69-
7069
Symbolics.derivative(::typeof(Ei), args::NTuple{1, Any}, ::Val{1}) = exp(args[1]) / args[1]
7170
Symbolics.derivative(::typeof(Si), args::NTuple{1, Any}, ::Val{1}) = sin(args[1]) / args[1]
7271
Symbolics.derivative(::typeof(Ci), args::NTuple{1, Any}, ::Val{1}) = cos(args[1]) / args[1]
7372
Symbolics.derivative(::typeof(Li), args::NTuple{1, Any}, ::Val{1}) = 1 / log(args[1])
74-
Symbolics.derivative(::typeof(Erfi), args::NTuple{1, Any}, ::Val{1}) = 2 / sqrt(2) * exp(args[1]^2)
73+
function Symbolics.derivative(::typeof(Erfi), args::NTuple{1, Any}, ::Val{1})
74+
2 / sqrt(2) * exp(args[1]^2)
75+
end
7576

7677
@syms 𝑥 si(𝑥) ci(𝑥) ei(𝑥) li(𝑥) erfi_(𝑥)
7778

@@ -82,9 +83,9 @@ guard_zero(x) = isequal(x, 0) ? one(x) : x
8283
function generate_homotopy(eq, x)
8384
eq = value(eq)
8485
x = value(x)
85-
86+
8687
if is_add(eq)
87-
return unique(union([generate_homotopy(t,x) for t in args(eq)]...))
88+
return unique(union([generate_homotopy(t, x) for t in args(eq)]...))
8889
end
8990

9091
p = transform(eq, x)
@@ -111,66 +112,69 @@ end
111112
@syms 𝛷(x)
112113

113114
partial_int_rules = [
114-
# trigonometric functions
115-
@rule 𝛷(sin(~x)) => (cos(~x) + si(~x), ~x)
116-
@rule 𝛷(cos(~x)) => (sin(~x) + ci(~x), ~x)
117-
@rule 𝛷(tan(~x)) => (log(cos(~x)), ~x)
118-
@rule 𝛷(csc(~x)) => (log(csc(~x) + cot(~x)) + log(sin(~x)), ~x)
119-
@rule 𝛷(sec(~x)) => (log(sec(~x) + tan(~x)) + log(cos(~x)), ~x)
120-
@rule 𝛷(cot(~x)) => (log(sin(~x)), ~x)
121-
# hyperbolic functions
122-
@rule 𝛷(sinh(~x)) => (cosh(~x), ~x)
123-
@rule 𝛷(cosh(~x)) => (sinh(~x), ~x)
124-
@rule 𝛷(tanh(~x)) => (log(cosh(~x)), ~x)
125-
@rule 𝛷(csch(~x)) => (log(tanh(~x / 2)), ~x)
126-
@rule 𝛷(sech(~x)) => (atan(sinh(~x)), ~x)
127-
@rule 𝛷(coth(~x)) => (log(sinh(~x)), ~x)
128-
# 1/trigonometric functions
129-
@rule 𝛷(1 / sin(~x)) => (log(csc(~x) + cot(~x)) + log(sin(~x)), ~x)
130-
@rule 𝛷(1 / cos(~x)) => (log(sec(~x) + tan(~x)) + log(cos(~x)), ~x)
131-
@rule 𝛷(1 / tan(~x)) => (log(sin(~x)) + log(tan(~x)), ~x)
132-
@rule 𝛷(1 / csc(~x)) => (cos(~x) + log(csc(~x)), ~x)
133-
@rule 𝛷(1 / sec(~x)) => (sin(~x) + log(sec(~x)), ~x)
134-
@rule 𝛷(1 / cot(~x)) => (log(cos(~x)) + log(cot(~x)), ~x)
135-
# 1/hyperbolic functions
136-
@rule 𝛷(1 / sinh(~x)) => (log(tanh(~x / 2)) + log(sinh(~x)), ~x)
137-
@rule 𝛷(1 / cosh(~x)) => (atan(sinh(~x)) + log(cosh(~x)), ~x)
138-
@rule 𝛷(1 / tanh(~x)) => (log(sinh(~x)) + log(tanh(~x)), ~x)
139-
@rule 𝛷(1 / csch(~x)) => (cosh(~x) + log(csch(~x)), ~x)
140-
@rule 𝛷(1 / sech(~x)) => (sinh(~x) + log(sech(~x)), ~x)
141-
@rule 𝛷(1 / coth(~x)) => (log(cosh(~x)) + log(coth(~x)), ~x)
142-
# inverse trigonometric functions
143-
@rule 𝛷(asin(~x)) => (~x * asin(~x) + sqrt(1 - ~x * ~x), ~x)
144-
@rule 𝛷(acos(~x)) => (~x * acos(~x) + sqrt(1 - ~x * ~x), ~x)
145-
@rule 𝛷(atan(~x)) => (~x * atan(~x) + log(~x * ~x + 1), ~x)
146-
@rule 𝛷(acsc(~x)) => (~x * acsc(~x) + atanh(1 - ^(~x, -2)), ~x)
147-
@rule 𝛷(asec(~x)) => (~x * asec(~x) + acosh(~x), ~x)
148-
@rule 𝛷(acot(~x)) => (~x * acot(~x) + log(~x * ~x + 1), ~x)
149-
# inverse hyperbolic functions
150-
@rule 𝛷(asinh(~x)) => (~x * asinh(~x) + sqrt(~x * ~x + 1), ~x)
151-
@rule 𝛷(acosh(~x)) => (~x * acosh(~x) + sqrt(~x * ~x - 1), ~x)
152-
@rule 𝛷(atanh(~x)) => (~x * atanh(~x) + log(~x + 1), ~x)
153-
@rule 𝛷(acsch(~x)) => (acsch(~x), ~x)
154-
@rule 𝛷(asech(~x)) => (asech(~x), ~x)
155-
@rule 𝛷(acoth(~x)) => (~x * acot(~x) + log(~x + 1), ~x)
156-
# logarithmic and exponential functions
157-
@rule 𝛷(log(~x)) => (~x + ~x * log(~x) + sum(pow_minus_rule(~x, -1); init = one(~x)), ~x)
158-
@rule 𝛷(1 / log(~x)) => (log(log(~x)) + li(~x), ~x)
159-
@rule 𝛷(exp(~x)) => (exp(~x) + ei(~x) + erfi_rule(~x), ~x)
160-
@rule 𝛷(^(exp(~x), ~k::is_neg)) => (^(exp(-~x), -~k), ~x)
161-
# square-root functions
162-
@rule 𝛷(^(~x, ~k::is_abs_half)) => (sum(sqrt_rule(~x, ~k); init = one(~x)), ~x);
163-
@rule 𝛷(sqrt(~x)) => (sum(sqrt_rule(~x, 0.5); init = one(~x)), ~x);
164-
@rule 𝛷(1 / sqrt(~x)) => (sum(sqrt_rule(~x, -0.5); init = one(~x)), ~x);
165-
# rational functions
166-
@rule 𝛷(1 / ^(~x::is_univar_poly, ~k::is_pos_int)) => (sum(pow_minus_rule(~x,-~k); init = one(~x)), ~x)
167-
@rule 𝛷(1 / ~x::is_univar_poly) => (sum(pow_minus_rule(~x, -1); init = one(~x)), ~x);
168-
@rule 𝛷(^(~x, -1)) => (log(~x), ~x)
169-
@rule 𝛷(^(~x, ~k::is_neg_int)) => (sum(^(~x, i) for i in (~k + 1):-1), ~x)
170-
@rule 𝛷(1 / ~x) => (log(~x), ~x)
171-
@rule 𝛷(^(~x, ~k::is_pos_int)) => (sum(^(~x, i + 1) for i in 1:(~k + 1)), ~x)
172-
@rule 𝛷(1) => (𝑥, 1)
173-
@rule 𝛷(~x) => ((~x + ^(~x, 2)), ~x)]
115+
# trigonometric functions
116+
@rule 𝛷(sin(~x)) => (cos(~x) + si(~x), ~x)
117+
@rule 𝛷(cos(~x)) => (sin(~x) + ci(~x), ~x)
118+
@rule 𝛷(tan(~x)) => (log(cos(~x)), ~x)
119+
@rule 𝛷(csc(~x)) => (log(csc(~x) + cot(~x)) + log(sin(~x)), ~x)
120+
@rule 𝛷(sec(~x)) => (log(sec(~x) + tan(~x)) + log(cos(~x)), ~x)
121+
@rule 𝛷(cot(~x)) => (log(sin(~x)), ~x)
122+
# hyperbolic functions
123+
@rule 𝛷(sinh(~x)) => (cosh(~x), ~x)
124+
@rule 𝛷(cosh(~x)) => (sinh(~x), ~x)
125+
@rule 𝛷(tanh(~x)) => (log(cosh(~x)), ~x)
126+
@rule 𝛷(csch(~x)) => (log(tanh(~x / 2)), ~x)
127+
@rule 𝛷(sech(~x)) => (atan(sinh(~x)), ~x)
128+
@rule 𝛷(coth(~x)) => (log(sinh(~x)), ~x)
129+
# 1/trigonometric functions
130+
@rule 𝛷(1 / sin(~x)) => (log(csc(~x) + cot(~x)) + log(sin(~x)), ~x)
131+
@rule 𝛷(1 / cos(~x)) => (log(sec(~x) + tan(~x)) + log(cos(~x)), ~x)
132+
@rule 𝛷(1 / tan(~x)) => (log(sin(~x)) + log(tan(~x)), ~x)
133+
@rule 𝛷(1 / csc(~x)) => (cos(~x) + log(csc(~x)), ~x)
134+
@rule 𝛷(1 / sec(~x)) => (sin(~x) + log(sec(~x)), ~x)
135+
@rule 𝛷(1 / cot(~x)) => (log(cos(~x)) + log(cot(~x)), ~x)
136+
# 1/hyperbolic functions
137+
@rule 𝛷(1 / sinh(~x)) => (log(tanh(~x / 2)) + log(sinh(~x)), ~x)
138+
@rule 𝛷(1 / cosh(~x)) => (atan(sinh(~x)) + log(cosh(~x)), ~x)
139+
@rule 𝛷(1 / tanh(~x)) => (log(sinh(~x)) + log(tanh(~x)), ~x)
140+
@rule 𝛷(1 / csch(~x)) => (cosh(~x) + log(csch(~x)), ~x)
141+
@rule 𝛷(1 / sech(~x)) => (sinh(~x) + log(sech(~x)), ~x)
142+
@rule 𝛷(1 / coth(~x)) => (log(cosh(~x)) + log(coth(~x)), ~x)
143+
# inverse trigonometric functions
144+
@rule 𝛷(asin(~x)) => (~x * asin(~x) + sqrt(1 - ~x * ~x), ~x)
145+
@rule 𝛷(acos(~x)) => (~x * acos(~x) + sqrt(1 - ~x * ~x), ~x)
146+
@rule 𝛷(atan(~x)) => (~x * atan(~x) + log(~x * ~x + 1), ~x)
147+
@rule 𝛷(acsc(~x)) => (~x * acsc(~x) + atanh(1 - ^(~x, -2)), ~x)
148+
@rule 𝛷(asec(~x)) => (~x * asec(~x) + acosh(~x), ~x)
149+
@rule 𝛷(acot(~x)) => (~x * acot(~x) + log(~x * ~x + 1), ~x)
150+
# inverse hyperbolic functions
151+
@rule 𝛷(asinh(~x)) => (~x * asinh(~x) + sqrt(~x * ~x + 1), ~x)
152+
@rule 𝛷(acosh(~x)) => (~x * acosh(~x) + sqrt(~x * ~x - 1), ~x)
153+
@rule 𝛷(atanh(~x)) => (~x * atanh(~x) + log(~x + 1), ~x)
154+
@rule 𝛷(acsch(~x)) => (acsch(~x), ~x)
155+
@rule 𝛷(asech(~x)) => (asech(~x), ~x)
156+
@rule 𝛷(acoth(~x)) => (~x * acot(~x) + log(~x + 1), ~x)
157+
# logarithmic and exponential functions
158+
@rule 𝛷(log(~x)) => (~x + ~x * log(~x) + sum(pow_minus_rule(~x, -1); init = one(~x)),
159+
~x);
160+
@rule 𝛷(1 / log(~x)) => (log(log(~x)) + li(~x), ~x)
161+
@rule 𝛷(exp(~x)) => (exp(~x) + ei(~x) + erfi_rule(~x), ~x)
162+
@rule 𝛷(^(exp(~x), ~k::is_neg)) => (^(exp(-~x), -~k), ~x)
163+
# square-root functions
164+
@rule 𝛷(^(~x, ~k::is_abs_half)) => (sum(sqrt_rule(~x, ~k); init = one(~x)), ~x);
165+
@rule 𝛷(sqrt(~x)) => (sum(sqrt_rule(~x, 0.5); init = one(~x)), ~x);
166+
@rule 𝛷(1 / sqrt(~x)) => (sum(sqrt_rule(~x, -0.5); init = one(~x)), ~x);
167+
# rational functions
168+
@rule 𝛷(1 / ^(~x::is_univar_poly, ~k::is_pos_int)) => (sum(pow_minus_rule(~x, -~k);
169+
init = one(~x)),
170+
~x);
171+
@rule 𝛷(1 / ~x::is_univar_poly) => (sum(pow_minus_rule(~x, -1); init = one(~x)), ~x);
172+
@rule 𝛷(^(~x, -1)) => (log(~x), ~x)
173+
@rule 𝛷(^(~x, ~k::is_neg_int)) => (sum(^(~x, i) for i in (~k + 1):-1), ~x)
174+
@rule 𝛷(1 / ~x) => (log(~x), ~x)
175+
@rule 𝛷(^(~x, ~k::is_pos_int)) => (sum(^(~x, i + 1) for i in 1:(~k + 1)), ~x)
176+
@rule 𝛷(1) => (𝑥, 1)
177+
@rule 𝛷(~x) => ((~x + ^(~x, 2)), ~x)]
174178

175179
function apply_partial_int_rules(eq, x)
176180
y, dy = Chain(partial_int_rules)(𝛷(value(eq)))
@@ -182,7 +186,7 @@ end
182186
function erfi_rule(eq)
183187
if is_univar_poly(eq)
184188
x = var(eq)
185-
return erfi_(x)
189+
return erfi_(x)
186190
end
187191
return 0
188192
end
@@ -227,13 +231,13 @@ end
227231

228232
function sqrt_rule(p, k)
229233
h = Any[p^k, p^(k + 1)]
230-
234+
231235
if !is_univar_poly(p)
232236
return h
233-
end
234-
237+
end
238+
235239
x = var(p)
236-
240+
237241
if poly_deg(p) == 2
238242
r, s = find_roots(p, x)
239243
l = leading(p, x)
@@ -255,4 +259,3 @@ function sqrt_rule(p, k)
255259
push!(h, log(0.5 * Δ + sqrt(p)))
256260
return h
257261
end
258-

0 commit comments

Comments
 (0)