Skip to content

Commit

Permalink
parametron file used for precompilation
Browse files Browse the repository at this point in the history
  • Loading branch information
jkosata committed Dec 13, 2022
1 parent 6e761d4 commit 67d9e96
Show file tree
Hide file tree
Showing 5 changed files with 21 additions and 28 deletions.
2 changes: 1 addition & 1 deletion src/DifferentialEquation.jl
Original file line number Diff line number Diff line change
Expand Up @@ -35,7 +35,7 @@ Return the dependent variables of `diff_eom`.
get_variables(diff_eom::DifferentialEquation) = collect(keys(diff_eom.equations))


is_harmonic(diff_eom::DifferentialEquation, t::Num) = all([is_harmonic(eq, t) for eq in values(diff_eom.equations)])
is_harmonic(diff_eom::DifferentialEquation, t::Num)::Bool = all([is_harmonic(eq, t) for eq in values(diff_eom.equations)])

"Pretty printing of the newly defined types"
function show_fields(object)
Expand Down
23 changes: 2 additions & 21 deletions src/HarmonicBalance.jl
Original file line number Diff line number Diff line change
Expand Up @@ -62,26 +62,7 @@ module HarmonicBalance
include("modules/LimitCycles.jl")
using .LimitCycles


@precompile_setup begin
# Putting some things in `setup` can reduce the size of the
# precompile file and potentially make loading faster.
#list = [OtherType("hello"), OtherType("world!")]

@variables Ω,γ,λ,F, x,θ,η,α, ω0, ω,t,T, ψ
@variables x(t)

natural_equation = d(d(x,t),t) + γ*d(x,t) + Ω^2*(1-λ*cos(2*ω*t+ψ))*x + α * x^3 +η *d(x,t) * x^2
forces = F*cos*t+θ)

@precompile_all_calls begin
# all calls in this block will be precompiled, regardless of whether
# they belong to your package or not (on Julia 1.8 and higher)
dEOM = DifferentialEquation(natural_equation + forces, x)
add_harmonic!(dEOM, x, ω)
harmonic_eq = get_harmonic_equations(dEOM);
end
end

precomp_path = (@__DIR__) * "/../test/parametron.jl"
@precompile_all_calls include(precomp_path)

end # module
2 changes: 1 addition & 1 deletion src/Symbolics_customised.jl
Original file line number Diff line number Diff line change
Expand Up @@ -95,7 +95,7 @@ exp_to_trig(x::Complex{Num}) = exp_to_trig(x.re)+im* exp_to_trig(x.im)

# sometimes, expressions get stored as Complex{Num} with no way to decode what real(x) and imag(x)
# this overloads the Num constructor to return a Num if x.re and x.im have similar arguments
Num(x::Complex{Num}) = isequal(x.re.val.arguments, x.im.val.arguments) ? Num(first(x.re.val.arguments)) : error("Cannot convert Complex{Num} " * string(x) * " to Num")
Num(x::Complex{Num})::Num = isequal(x.re.val.arguments, x.im.val.arguments) ? Num(first(x.re.val.arguments)) : error("Cannot convert Complex{Num} " * string(x) * " to Num")


#=
Expand Down
12 changes: 12 additions & 0 deletions test/firstplot_benchmark.jl
Original file line number Diff line number Diff line change
@@ -0,0 +1,12 @@
using HarmonicBalance

@variables α, ω, ω0, F, t, η, x(t) # declare constant variables and a function x(t)
diff_eq = DifferentialEquation(d(x,t,2) + ω0^2*x + α*x^3 + η*d(x,t)*x^2 ~ F*cos*t), x)
add_harmonic!(diff_eq, x, ω) # specify the ansatz x = u(T) cos(ωt) + v(T) sin(ωt)

# implement ansatz to get harmonic equations
harmonic_eq = get_harmonic_equations(diff_eq)

fixed ==> 1., ω0 => 1.0, F => 0.01, η=>0.1) # fixed parameters
varied = ω => LinRange(0.9, 1.2, 100) # range of parameter values
result = get_steady_states(harmonic_eq, varied, fixed, show_progress=false);
10 changes: 5 additions & 5 deletions test/parametron.jl
Original file line number Diff line number Diff line change
@@ -1,6 +1,6 @@
using HarmonicBalance
using Symbolics
using Test
#using Test # do not use Test as this file is used for precompilation

@variables Ω,γ,λ,F, x,θ,η,α, ω0, ω,t,T, ψ
@variables x(t)
Expand All @@ -15,7 +15,7 @@ p = HarmonicBalance.Problem(harmonic_eq);

fixed ==> 1.0=> 1E-2, λ => 5E-2, F => 1E-3, α => 1., η=>0.3, θ => 0, ψ => 0)
varied = ω => LinRange(0.9, 1.1, 100)
res = get_steady_states(p, varied, fixed)
res = get_steady_states(p, varied, fixed, show_progress=false)

# save the result, try and load in the next step
#current_path = @__DIR__
Expand All @@ -24,7 +24,7 @@ res = get_steady_states(p, varied, fixed)
# try to run a 2D calculation
fixed ==> 1.0=> 1E-2, F => 1E-3, α => 1., η=>0.3, θ => 0, ψ => 0)
varied ==> LinRange(0.9, 1.1, 10), λ => LinRange(0.01, 0.05, 10))
res = get_steady_states(p, varied, fixed)
res = get_steady_states(p, varied, fixed, show_progress=false)


###
Expand All @@ -35,5 +35,5 @@ ref1 = (Ω^2)*u1 + F*cos(θ) + γ*Differential(T)(u1) + (3//4)*α*(u1^3) + γ*ω
ref2 = γ*Differential(T)(v1) +^2)*v1 + (3//4)*α*(v1^3) + (3//4)*α*(u1^2)*v1 + (1//4)*η*(u1^2)*Differential(T)(v1) + (3//4)*η*(v1^2)*Differential(T)(v1) + (1//2)*λ*^2)*v1*cos(ψ) + (1//2)*η*u1*v1*Differential(T)(u1) + (1//2)*λ*^2)*u1*sin(ψ) - F*sin(θ) -^2)*v1 - (2//1)*ω*Differential(T)(u1) - γ*ω*u1 - (1//4)*η*ω*(u1^3) - (1//4)*η*ω*(v1^2)*u1

averaged = HarmonicBalance._remove_brackets(harmonic_eq)
@test isequal(simplify(expand(averaged[1]-ref1)), 0)
@test isequal(simplify(expand(averaged[2]-ref2)), 0)
@assert isequal(simplify(expand(averaged[1]-ref1)), 0)
@assert isequal(simplify(expand(averaged[2]-ref2)), 0)

0 comments on commit 67d9e96

Please sign in to comment.