Skip to content

Use Symbolics #825

New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Merged
merged 8 commits into from
Mar 1, 2021
Merged
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
6 changes: 4 additions & 2 deletions Project.toml
Original file line number Diff line number Diff line change
@@ -1,7 +1,7 @@
name = "ModelingToolkit"
uuid = "961ee093-0014-501f-94e3-6117800e7a78"
authors = ["Chris Rackauckas <accounts@chrisrackauckas.com>"]
version = "5.9.1"
version = "5.10.0"

[deps]
ArrayInterface = "4fba245c-0d91-5ea0-9b3e-6abc04ee57a9"
Expand Down Expand Up @@ -33,6 +33,7 @@ SparseArrays = "2f01184e-e22b-5df5-ae63-d93ebab69eaf"
SpecialFunctions = "276daf66-3868-5448-9aa4-cd146d93841b"
StaticArrays = "90137ffa-7385-5640-81b9-e52037218182"
SymbolicUtils = "d1185830-fcd6-423d-90d6-eec64667417b"
Symbolics = "0c5d862f-8b57-4792-8d23-62f2024744c7"
TreeViews = "a2a6695c-b41b-5b7d-aed9-dbfdeacea5d7"
UnPack = "3a884ed6-31ef-47d7-9d2a-63182c4928ed"
Unitful = "1986cc42-f94f-5a68-af5c-568840ba703d"
Expand Down Expand Up @@ -62,7 +63,8 @@ SciMLBase = "1.3"
Setfield = "0.7"
SpecialFunctions = "0.7, 0.8, 0.9, 0.10, 1.0"
StaticArrays = "0.10, 0.11, 0.12, 1.0"
SymbolicUtils = "0.8"
SymbolicUtils = "0.8.3"
Symbolics = "0.1"
TreeViews = "0.3"
UnPack = "0.1, 1.0"
Unitful = "1.1"
Expand Down
175 changes: 25 additions & 150 deletions src/ModelingToolkit.jl
Original file line number Diff line number Diff line change
Expand Up @@ -25,159 +25,45 @@ RuntimeGeneratedFunctions.init(@__MODULE__)
using RecursiveArrayTools

import SymbolicUtils
import SymbolicUtils: Term, Add, Mul, Pow, Sym, FnType,
@rule, Rewriters, substitute, similarterm,
promote_symtype

import SymbolicUtils: istree, arguments, operation, similarterm, promote_symtype,
Symbolic, Term, Add, Mul, Pow, Sym, FnType,
@rule, Rewriters, substitute
using SymbolicUtils.Code
import SymbolicUtils.Code: toexpr

import SymbolicUtils.Rewriters: Chain, Postwalk, Prewalk, Fixpoint

using Reexport
@reexport using Symbolics
export @derivatives
using Symbolics: _parse_vars, value, makesym, @derivatives, get_variables,
exprs_occur_in, solve_for, build_expr
import Symbolics: rename, get_variables!, _solve, hessian_sparsity,
jacobian_sparsity, islinear, _iszero, _isone,
tosymbol, lower_varname, diff2term, var_from_nested_derivative,
BuildTargets, JuliaTarget, StanTarget, CTarget, MATLABTarget,
ParallelForm, SerialForm, MultithreadedForm, build_function,
unflatten_long_ops

import DiffEqBase: @add_kwonly
using LinearAlgebra: LU, BlasInt

import LightGraphs: SimpleDiGraph, add_edge!

import TreeViews

using Requires

export Num, Variable
"""
$(TYPEDEF)

Wrap anything in a type that is a subtype of Real
"""
struct Num <: Real
val
end

const show_numwrap = Ref(false)

Num(x::Num) = x # ideally this should never be called
(n::Num)(args...) = Num(value(n)(map(value,args)...))
value(x) = x
value(x::Num) = x.val

SciMLBase.issymbollike(::Num) = true
SciMLBase.issymbollike(::SymbolicUtils.Symbolic) = true

SymbolicUtils.@number_methods(Num,
Num(f(value(a))),
Num(f(value(a), value(b))))

for C in [Complex, Complex{Bool}]
for fun in [:toexpr]
@eval begin
Base.:*(x::Num, z::$C) = Complex(x * real(z), x * imag(z))
Base.:*(z::$C, x::Num) = Complex(real(z) * x, imag(z) * x)
end
end

Base.:+(x::Num, z::Complex) = Complex(x + real(z), imag(z))
Base.:+(z::Complex, x::Num) = Complex(real(z) + x, imag(z))
Base.:-(x::Num, z::Complex) = Complex(x - real(z), -imag(z))
Base.:-(z::Complex, x::Num) = Complex(real(z) - x, imag(z))

function Base.inv(z::Complex{Num})
a, b = reim(z)
den = a^2 + b^2
Complex(a/den, -b/den)
end
function Base.:/(x::Complex{Num}, y::Complex{Num})
a, b = reim(x)
c, d = reim(y)
den = c^2 + d^2
Complex((a*c + b*d)/den, (b*c - a*d)/den)
end
function $fun(eq::Equation; kw...)
Expr(:(=), $fun(eq.lhs; kw...), $fun(eq.rhs; kw...))
end

function Base.show(io::IO, z::Complex{<:Num})
r, i = reim(z)
compact = get(io, :compact, false)
show(io, r)
print(io, (compact ? "+" : " + ") * "(")
show(io, i)
print(io, ")*im")
end

SymbolicUtils.simplify(n::Num; kw...) = Num(SymbolicUtils.simplify(value(n); kw...))

SymbolicUtils.symtype(n::Num) = symtype(n.val)

function Base.iszero(x::Num)
_x = SymbolicUtils.to_mpoly(value(x))[1]
return (_x isa Number || _x isa SymbolicUtils.MPoly) && iszero(_x)
end

import SymbolicUtils: <ₑ, Symbolic, Term, operation, arguments

Base.show(io::IO, n::Num) = show_numwrap[] ? print(io, :(Num($(value(n))))) : Base.show(io, value(n))

Base.promote_rule(::Type{<:Number}, ::Type{<:Num}) = Num
Base.promote_rule(::Type{<:Symbolic{<:Number}}, ::Type{<:Num}) = Num
function Base.getproperty(t::Union{Add, Mul, Pow, Term}, f::Symbol)
if f === :op
Base.depwarn("`x.op` is deprecated, use `operation(x)` instead", :getproperty, force=true)
operation(t)
elseif f === :args
Base.depwarn("`x.args` is deprecated, use `arguments(x)` instead", :getproperty, force=true)
arguments(t)
else
getfield(t, f)
end
end
<ₑ(s::Num, x) = value(s) <ₑ value(x)
<ₑ(s, x::Num) = value(s) <ₑ value(x)
<ₑ(s::Num, x::Num) = value(s) <ₑ value(x)

for T in (Integer, Rational)
@eval Base.:(^)(n::Num, i::$T) = Num(value(n)^i)
end

macro num_method(f, expr, Ts=nothing)
if Ts === nothing
Ts = [Any]
else
@assert Ts.head == :tuple
# e.g. a tuple or vector
Ts = Ts.args
$fun(eqs::AbstractArray; kw...) = map(eq->$fun(eq; kw...), eqs)
$fun(x::Integer; kw...) = x
$fun(x::AbstractFloat; kw...) = x
end

ms = [quote
$f(a::$T, b::$Num) = $expr
$f(a::$Num, b::$T) = $expr
end for T in Ts]
quote
$f(a::$Num, b::$Num) = $expr
$(ms...)
end |> esc
end

"""
tosymbolic(a::Union{Sym,Num}) -> Sym{Real}
tosymbolic(a::T) -> T
"""
tosymbolic(a::Num) = tosymbolic(value(a))
tosymbolic(a::Sym) = tovar(a)
tosymbolic(a) = a
@num_method Base.isless (val = isless(tosymbolic(a), tosymbolic(b)); val isa Bool ? val : Num(val)) (Real,)
@num_method Base.:(<) (val = tosymbolic(a) < tosymbolic(b) ; val isa Bool ? val : Num(val)) (Real,)
@num_method Base.:(<=) (val = tosymbolic(a) <= tosymbolic(b) ; val isa Bool ? val : Num(val)) (Real,)
@num_method Base.:(>) (val = tosymbolic(a) > tosymbolic(b) ; val isa Bool ? val : Num(val)) (Real,)
@num_method Base.:(>=) (val = tosymbolic(a) >= tosymbolic(b) ; val isa Bool ? val : Num(val)) (Real,)
@num_method Base.:(==) (val = tosymbolic(a) == tosymbolic(b) ; val isa Bool ? val : Num(val)) (AbstractFloat,Number)
@num_method Base.isequal isequal(tosymbolic(a), tosymbolic(b)) (AbstractFloat, Number, Symbolic)

Base.hash(x::Num, h::UInt) = hash(value(x), h)

Base.convert(::Type{Num}, x::Symbolic{<:Number}) = Num(x)
Base.convert(::Type{Num}, x::Number) = Num(x)
Base.convert(::Type{Num}, x::Num) = x

Base.convert(::Type{<:Array{Num}}, x::AbstractArray) = map(Num, x)
Base.convert(::Type{<:Array{Num}}, x::AbstractArray{Num}) = x
Base.convert(::Type{Sym}, x::Num) = value(x) isa Sym ? value(x) : error("cannot convert $x to Sym")

LinearAlgebra.lu(x::Array{Num}; check=true, kw...) = sym_lu(x; check=check)

"""
$(TYPEDEF)

Expand Down Expand Up @@ -211,16 +97,10 @@ include("bipartite_graph.jl")
using .BipartiteGraphs

include("variables.jl")
include("context_dsl.jl")
include("differentials.jl")
include("parameters.jl")

include("equations.jl")
include("utils.jl")
include("linearity.jl")
include("solve.jl")
include("direct.jl")
include("domains.jl")
include("register_function.jl")

include("systems/abstractsystem.jl")

Expand Down Expand Up @@ -252,10 +132,6 @@ include("systems/alias_elimination.jl")
include("structural_transformation/StructuralTransformations.jl")
@reexport using .StructuralTransformations

include("latexify_recipes.jl")
include("build_function.jl")
include("extra_functions.jl")

export ODESystem, ODEFunction, ODEFunctionExpr, ODEProblemExpr
export SDESystem, SDEFunction, SDEFunctionExpr, SDESystemExpr
export SystemStructure
Expand Down Expand Up @@ -296,7 +172,6 @@ export asgraph, asdigraph
export toexpr, get_variables
export simplify, substitute
export build_function
export @register
export modelingtoolkitize
export @variables, @parameters
export @named
Expand Down
Loading