forked from JuliaSymbolics/SymbolicIntegration.jl
-
Notifications
You must be signed in to change notification settings - Fork 0
Expand file tree
/
Copy pathfrontend.jl
More file actions
117 lines (104 loc) · 3.99 KB
/
Copy pathfrontend.jl
File metadata and controls
117 lines (104 loc) · 3.99 KB
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
include("string_manipulation_helpers.jl")
"""
Applies iteratively rules from the RULES array until a result is found.
returns a tuple:
if found a rule to apply, (solution, true)
if not, (original problem, false)
"""
function apply_rule(problem)
result = nothing
for (i, rule) in enumerate(RULES)
result = rule(problem)
if result !== nothing
if result===problem
VERBOSE && println("Infinite cycle created by rule $(IDENTIFIERS[i]) applied on ", problem)
continue
end
if VERBOSE && !in(IDENTIFIERS[i], SILENCE)
s = pretty_print_rule(rule, IDENTIFIERS[i])
printstyled("┌-------Applied rule $(IDENTIFIERS[i]) on ";);
printstyled(problem; color = :light_red)
for ss in split(s, '\n')
printstyled("\n| ";); printstyled(ss;bold=true)
end
printstyled("\n└-------with result: ";)
printstyled(result, "\n"; color = :light_blue)
end
in(IDENTIFIERS[i], SILENCE) && pop!(SILENCE)
return (result, true)
end
end
VERBOSE && println("No rule found for ", problem)
return (problem, false)
end
"""
ins = inverse not simplified
This function creates a term with negative power that doesnt simplify
automatically to a division, like would happen with the ^ function
```
julia> SymbolicIntegration.ins(Symbolics.unwrap(x))
x^-1
julia> SymbolicIntegration.ins(Symbolics.unwrap(x^3))
x^-3
```
"""
function ins(expr)
t = (@rule (~u)^(~m) => ~)(expr)
t!==nothing && return SymbolicUtils.Term{Number}(^,[t[:u],-t[:m]])
return SymbolicUtils.Term{Number}(^,[expr,-1])
end
# TODO add threaded for speed?
function repeated_prewalk(expr)
!iscall(expr) && return expr
if operation(expr)===∫
(new_expr,success) = apply_rule(expr)
# r1 and r2 are needed bc of neim problem
if !success
r2 = @rule ∫((~n)/*(~~d),~x) => ∫(~n*prod([ins(el) for el in ~~d]),~x)
r2r = r2(expr)
if r2r!==nothing
VERBOSE && println("integration of ", expr, " failed, trying with this mathematically equivalent integrand:\n$r2r")
(new_expr,success) = apply_rule(r2r)
if success && new_expr===expr
success=false
end
end
end
if !success
r1 = @rule ∫((~n)/(~d),~x) => ∫(~n*ins(~d),~x)
r1r = r1(expr)
if r1r!==nothing
VERBOSE && println("integration of ", expr, " failed, trying with this mathematically equivalent integrand:\n$r1r")
(new_expr,success) = apply_rule(r1r)
# if success we know r1r!=new_expr
# but clud be new_expr==expr
if success && new_expr===expr
success=false
end
end
end
if !success
# TODO Can this be a bad idea sometimes?
simplified_expr = simplify(expr, expand=true)
if simplified_expr !== expr
VERBOSE && println("integration of ", expr, " failed, trying with the expanded version:\n", simplified_expr)
(new_expr,success) = apply_rule(simplified_expr)
end
end
success && return repeated_prewalk(new_expr)
end
expr = SymbolicUtils.maketerm(
typeof(expr),
operation(expr),
map(repeated_prewalk, arguments(expr)),
SymbolicUtils.metadata(expr)
)
return expr
end
function integrate_rule_based(integrand::Symbolics.Num, int_var::Symbolics.Num; use_gamma::Bool=false, verbose::Bool=true, kwargs...)
global VERBOSE
VERBOSE = verbose
return simplify(repeated_prewalk(∫(integrand,int_var)))
end
integrate_rule_based(integrand::SymbolicUtils.BasicSymbolic{Real}, int_var::SymbolicUtils.BasicSymbolic{Real}; kwargs...) =
integrate_rule_based(Num(integrand), Num(int_var); kwargs...)