JuliaDynamics/DynamicalSystems.jl

produce_orbitdiagram problem.

gacarita opened this issue · 7 comments

Hi there.

I'm trying some features in the DynamicalSystems, however i'm facing some issues.

I'm trying to use produce_orbitdiagram function. Everything sounds good such defining the dynamicalsystem, evolving the trajectory and constructing the poincare section. When i try to construct the orbitdiagram i'm facing issues.


.
.
.

ds = ContinuousDynamicalSystem(eom!, u0, p)

plane = (3, 0.0)
pvalues = collect(0.1:0.15:201)
i = Int(1)
plane = (3, 0.0)
tf = 200.0
p_index = 1

output = produce_orbitdiagram(ds, plane, i,pvalues;
                              tfinal = tf, Ttr = 200.0)


I'm getting this error :

MethodError: no method matching produce_orbitdiagram(::ContinuousDynamicalSystem{true, Vector{Float64}, 4, typeof(eom!), Float64, DynamicalSystemsBase.var"#6#12"{ForwardDiff.JacobianConfig{ForwardDiff.Tag{DynamicalSystemsBase.var"#5#11"{typeof(eom!), Float64, Float64}, Float64}, Float64, 4, Tuple{Vector{ForwardDiff.Dual{ForwardDiff.Tag{DynamicalSystemsBase.var"#5#11"{typeof(eom!), Float64, Float64}, Float64}, Float64, 4}}, Vector{ForwardDiff.Dual{ForwardDiff.Tag{DynamicalSystemsBase.var"#5#11"{typeof(eom!), Float64, Float64}, Float64}, Float64, 4}}}}, DynamicalSystemsBase.var"#5#11"{typeof(eom!), Float64, Float64}, Vector{Float64}}, Matrix{Float64}, true}, ::Tuple{Int64, Float64}, ::Int64, ::Vector{Float64}; tfinal=200.0, Ttr=200.0)
Closest candidates are:
  produce_orbitdiagram(::ContinuousDynamicalSystem{IIP, S, D, F, P, JAC, JM, IAD} where {F, P, JAC, JM, IAD}, ::Any, ::Any, ::Any, !Matched::Any; tfinal, direction, printparams, warning, Ttr, u0, rootkw, diffeq...) where {IIP, S, D} at C:\Users\gabri\.julia\packages\ChaosTools\OwxTS\src\orbitdiagrams\poincare.jl:318
in eval at [base\boot.jl:360](https://github.com/JuliaDynamics/DynamicalSystems.jl/issues/new?assignees=&labels=&template=bug_report.md&title=#) 
in top-level scope at [cr3bp_chaostools.jl:96](https://github.com/JuliaDynamics/DynamicalSystems.jl/issues/new?assignees=&labels=&template=bug_report.md&title=#)

I cannot understand where is the error. Any help ?

Gabriel.

Hi,
Please provide a Minimal Working Example, i.e., a code that can run and produce this error. This helps me produce the error.

(also, please do not ignore the issue template. It is there for a reason: to solve your bug faster. This means, you should also report the versions that you used)

Here is the code i'm using

using DataFrames, DifferentialEquations, DynamicalSystems, Plots

@inline @inbounds function vy0(C,p,v)
	"""
	Initial Velocity for cte jacobi conservation
		in:
			C = cte jacobi
			v = position and velocities of the particle
			cartesian coordinates
			p = mass ratio parameter == μ2
		out:
			ydot = velocity in y

	"""
		x, y, x̂= v
		μ = p
		#distance vectors
		r1 = sqrt((x+μ)^2 + y^2)
		r2 = sqrt((x - (1-μ))^2 + y^2)
		ŷ = sqrt(x^2 + y^2 - x̂^2 + 2*((1-μ)/r1 + μ/r2) - C)
		if x > 0
			return -ŷ
		end
		if x < 0
			return ŷ
		end
end


@inline @inbounds function eom!(du,u,p,t)
    """
    Planar Circular Restricted 3body problem
    Equations of motion in co-rotating frame

        in:
            u = position and velocities of the particle
            cartesian coordinates
                x = u[1]
                y = u[2]
                vx = u[3]
                vy = u[4]
            p = mass ratio parameter == μ2
            t = time span
        out:
            du = integration of the acel. and vel.

    """
        local μ = p

        #cartesian coordinates
        x = u[1]
        y = u[2]
        vx = u[3]
        vy = u[4]

        r1 = sqrt((x + μ)^2 + y^2)
        r2 = sqrt((x - (1-μ))^2 + y^2)

        du[1] = vx # F1
        du[2] = vy # F2
        du[3] =  2*vy  + x - (1 - μ)*(μ + x)/(r1^3) - μ*(x -(1-μ))/(r2^3) # F3
        du[4] = -2*vx  + y - ((1 - μ)/(r1^3) + μ/(r2^3))*y # F4
        return
    end

p=0.07

x0=-1.24
y0=0.19
vx0=0
vyu = vy0(-1.2,p,[x0,y0,vx0])

u0 = [x0,y0,vx0,vyu]
plane = (3, 0.0)

ds = ContinuousDynamicalSystem(eom!, u0, p)
psos = poincaresos(ds, plane, 6000.0; u0 = u0,direction = 0,rootkw = (xrtol = 1e-9, atol = 1e-9),
        diffeq = (alg = Vern9()))
scatter(psos[:, 1], psos[:, 2]; markersize = 2.0)

#code works fine until here

pvalues = collect(0.1:0.15:201)
i = Int(1)
plane = (3, 0.0)
tf = 200.0
p_index = 1

output = produce_orbitdiagram(ds, plane, i,pvalues;
                              tfinal = tf, Ttr = 200.0)


The call signature for produce_orbitdiagram is:
image

You only provide 4 positional arguments, you miss p_index.

Sorry i missed that part when i shared the code, however even setting p_index=1 in to the output

output = produce_orbitdiagram(ds, plane, i,p_index,pvalues;
                              tfinal = tf, Ttr = 200.0)

I still get an error

no method matching setindex!(::Float64, ::Float64, ::Int64)

I see. Well it's because your parameter "container" isn't a container, but a number, p=0.07. Do p = [0.07] instead. Numbers are not mutable, but vectors are.

Gotcha. Thank you.