Marcel Rüdiger
Marcel Rüdiger

Reputation: 13

How to Plot a function of two variables in Julia with pyplot

I'm trying to plot a function of two variables with pyplot in Julia. The working starting-point is the following (found here at StackOverflow):

function f(z,t)
    return z*t
end

z = linspace(0,5,11)
t = linspace(0,40,4)

for tval in t
   plot(z, f(z, tval))
end
show()

This works right for me and is giving me exactly what I wanted: a field of lines.

My own functions are as follows:

## needed functions ##

const gamma_0 = 6
const Ksch = 1.2
const Kver = 1.5

function Kvc(vc)
    if vc <= 0 
        return 0
    elseif vc < 20 
        return (100/vc)^0.1
    elseif vc < 100
        return 2.023/(vc^0.153)
    elseif vc == 100
        return 1
    elseif vc > 100
        return 1.380/(vc^0.07)
    else 
        return 0
    end
end 

function Kgamma(gamma_t)
    return 1-((gamma_t-gamma_0)/100)
end

function K(gamma_t, vc)
    return Kvc(vc)*Kgamma(gamma_t)*Ksch*Kver
end

I've tried to plot them as follows:

i = linspace(0,45,10)
j = linspace(0,200,10)
for i_val in i
    plot(i,K(i,j))
end

This gives me the following Error:

isless has no method matching isless(::Int64, ::Array{Float64,1}) while loading In[51], in expression starting on line 3

in Kvc at In[17]:2 in anonymous at no file:4

Obviously, my function cant deal with an array.

Next try:

i = linspace(0,200,11)
j = linspace(0,45,11)
for i_val in i
    plot(i_val,map(K,i_val,j))
end

gives me a empty plot only with axes

Can anybody please give me a hint...

EDIT

A simpler example:

using PyPlot
function P(n,M)
    return (M*n^3)/9550
end
M = linspace(1,5,5)
n = linspace(0,3000,3001)

for M_val in M
    plot(n,P(n,M_val))
end
show()

Solution

OK, with your help I found this solution for the shortened example which works for me as intended:

function P(n,M)
    result = Array(Float64, length(n))
    for (idx, val) in enumerate(n)
        result[idx] = (M*val^3)/9550
    end
    return result
end
n = linspace(0,3000,3001)

for M_val = 1:5
    plot(n,P(n,M_val))
end
show()

This gives me what I wanted for this shortened example. The remainig question is: could it be done in a simpler more elegant way?

I'll try to apply it to the original example and post it when I'll succed.

Upvotes: 0

Views: 2424

Answers (2)

Chris Rackauckas
Chris Rackauckas

Reputation: 19162

< is defined for scalars. I think you need to broadcast it for arrays, i.e. use .<. Example:

julia> x = 2
2

julia> x < 3
true

julia> x < [3 4]
ERROR: MethodError: no method matching isless(::Int64, ::Array{Int64,2})
Closest candidates are:
  isless(::Real, ::AbstractFloat)
  isless(::Real, ::Real)
  isless(::Integer, ::Char)
 in <(::Int64, ::Array{Int64,2}) at .\operators.jl:54
 in eval(::Module, ::Any) at .\boot.jl:234
 in macro expansion at .\REPL.jl:92 [inlined]
 in (::Base.REPL.##1#2{Base.REPL.REPLBackend})() at .\event.jl:46

julia> x .< [3 4]
1x2 BitArray{2}:
 true  true

Upvotes: 1

Michael Ohlrogge
Michael Ohlrogge

Reputation: 10990

I don't completely follow all the details of what you're trying to accomplish, but here are examples on how you can modify a couple of your functions so that they accept and return arrays:

function Kvc(vc)
    result = Array(Float64, length(vc))
    for (idx, val) in enumerate(vc)
        if val <= 0
            result[idx] = 0
        elseif val < 20
            result[idx] = (100/val)^0.1
        elseif val < 100
            result[idx] = 2.023/(val^0.153)
        elseif val == 100
            result[idx] = 1
        elseif val > 100
            result[idx] = 1.380/(val^0.07)
        else
            result[idx] = 0
        end
    end
    return result
end 

function Kgamma(gamma_t)
    return ones(length(gamma_t))-((gamma_t - gamma_0)/100)
end

Also, for your loop, I think you probably want something like:

for i_val in i
    plot(i_val,K(i_val,j))
end

rather than plot(i, K(i,j), as that would just print the same thing over and over.

Upvotes: 1

Related Questions