notebook.community

Edit and run

Exercise1

線形補間をするコードを書いてみる. 詳細: https://github.com/myuuuuun/oyama_seminar2016/tree/master/exercise/ex01





In [1]:



    


# imports
using Gadfly
set_default_plot_size(24cm, 12cm)



    


















    






WARNING: Method definition describe(AbstractArray{T, N} where N where T) in module StatsBase at /Users/akira/.julia/v0.7/StatsBase/src/scalarstats.jl:573 overwritten in module DataFrames at /Users/akira/.julia/v0.7/DataFrames/src/abstractdataframe/abstractdataframe.jl:407.

searchsorted() の使い方





In [2]:



    


array = collect(0:2:10)
println(array)

# Return the index of the first value in a greater than or equal to x
println("first, 3.5: ", searchsortedfirst(array, 3.5))
println("first, -2: ", searchsortedfirst(array, -2))
println("first, 11: ", searchsortedfirst(array, 11))

# Return the index of the first value in a less than or equal to x
println("last, 3.5: ", searchsortedlast(array, 3.5))
println("last, -2: ", searchsortedlast(array, -2))
println("last, 11: ", searchsortedlast(array, 11))

# Returns the range of indices of a which compare as equal to x
println("range, 3.5: ", searchsorted(array, 3.5))



    


















    






[0, 2, 4, 6, 8, 10]
first, 3.5: 3
first, -2: 1
first, 11: 7
last, 3.5: 2
last, -2: 0
last, 11: 6
range, 3.5:

function in function による実装





In [3]:



    


function my_lin_interp(grid, vals)
    function func(x)
        lower_index = searchsortedlast(grid, x)
        upper_index = lower_index + 1
        if lower_index == 0 || lower_index == length(grid)
            throw(DomainError())
        end

        grid_step = grid[upper_index] - grid[lower_index]
        val_step = vals[upper_index] - vals[lower_index]
        iterp_val = vals[lower_index] + (val_step / grid_step)* (x - grid[lower_index])
        
        return iterp_val
    end

    return func
end



    


















    
Out[3]:








my_lin_interp (generic function with 1 method)

















In [4]:



    


grid = collect(0:2:10)
val = collect(0:2:10)
interp_func = my_lin_interp(grid, val)
println(interp_func(0))
println(interp_func(1))
println(interp_func(2))



    


















    






0.0
1.0
2.0

$y = cos(x), -3\pi \leq x \leq 3\pi$ を補間してみる





In [5]:



    


grid2 = linspace(-3pi, 3pi, 201)
val2 = cos.(grid2)

grid3 = linspace(-3pi, 3pi, 19)
val3 = cos.(grid3)

set_default_plot_size(24cm, 6cm)

plot(
    layer(x=grid2, y=val2, Geom.line, Theme(line_width=2pt)),
    layer(x=grid3, y=val3, Geom.point, Geom.line, Theme(default_color=colorant"red"), order=2),
    Guide.xlabel("x"), Guide.ylabel("cosx"), Guide.title("Interpolation of cos x")
    )



    


















    
Out[5]:









  
    x
  
  
    -35
    -30
    -25
    -20
    -15
    -10
    -5
    0
    5
    10
    15
    20
    25
    30
    35
    -30
    -29
    -28
    -27
    -26
    -25
    -24
    -23
    -22
    -21
    -20
    -19
    -18
    -17
    -16
    -15
    -14
    -13
    -12
    -11
    -10
    -9
    -8
    -7
    -6
    -5
    -4
    -3
    -2
    -1
    0
    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
    -40
    -20
    0
    20
    40
    -30
    -28
    -26
    -24
    -22
    -20
    -18
    -16
    -14
    -12
    -10
    -8
    -6
    -4
    -2
    0
    2
    4
    6
    8
    10
    12
    14
    16
    18
    20
    22
    24
    26
    28
    30
  

  
    
      
    
    
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
    
    
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
      
    
    
      
        
      
      
        
      
      
        
          
          
          
          
          
          
          
          
          
          
          
          
          
          
          
          
          
          
          
        
      
    
    
      
        
        
      
      
        
      
      
        
      
      
        
        
      
    
  

  
    -3.5
    -3.0
    -2.5
    -2.0
    -1.5
    -1.0
    -0.5
    0.0
    0.5
    1.0
    1.5
    2.0
    2.5
    3.0
    3.5
    -3.0
    -2.9
    -2.8
    -2.7
    -2.6
    -2.5
    -2.4
    -2.3
    -2.2
    -2.1
    -2.0
    -1.9
    -1.8
    -1.7
    -1.6
    -1.5
    -1.4
    -1.3
    -1.2
    -1.1
    -1.0
    -0.9
    -0.8
    -0.7
    -0.6
    -0.5
    -0.4
    -0.3
    -0.2
    -0.1
    0.0
    0.1
    0.2
    0.3
    0.4
    0.5
    0.6
    0.7
    0.8
    0.9
    1.0
    1.1
    1.2
    1.3
    1.4
    1.5
    1.6
    1.7
    1.8
    1.9
    2.0
    2.1
    2.2
    2.3
    2.4
    2.5
    2.6
    2.7
    2.8
    2.9
    3.0
    -4
    -2
    0
    2
    4
    -3.0
    -2.8
    -2.6
    -2.4
    -2.2
    -2.0
    -1.8
    -1.6
    -1.4
    -1.2
    -1.0
    -0.8
    -0.6
    -0.4
    -0.2
    0.0
    0.2
    0.4
    0.6
    0.8
    1.0
    1.2
    1.4
    1.6
    1.8
    2.0
    2.2
    2.4
    2.6
    2.8
    3.0
  
  
    cosx
  
  
    Interpolation of cos x
  


  
  

  
    
  





















In [6]:



    


interp_func2 = my_lin_interp(grid3, val3)
@printf "cos(0) の 真値: %.5f, 補間値: %.5f\n" cos(0) interp_func2(0)
@printf "cos(π/4) の 真値: %.5f, 補間値: %.5f\n" cos(pi/4) interp_func2(pi/4)
@printf "cos(π/3) の 真値: %.5f, 補間値: %.5f\n" cos(pi/3) interp_func2(pi/3)

println("\n参考: 1/√2 = ", 1/1.4142)



    


















    






cos(0) の 真値: 1.00000, 補間値: 1.00000
cos(π/4) の 真値: 0.70711, 補間値: 0.62500
cos(π/3) の 真値: 0.50000, 補間値: 0.50000

参考: 1/√2 = 0.7071135624381276

type による実装





In [7]:



    


immutable MyInterp
    grid::Array
    vals::Array
end

function (points::MyInterp)(x)
    grid = points.grid; vals = points.vals
    lower_index = searchsortedlast(points.grid, x)
    upper_index = lower_index + 1
    if lower_index == 0 || lower_index == length(grid)
        throw(DomainError())
    end

    grid_step = grid[upper_index] - grid[lower_index]
    val_step = vals[upper_index] - vals[lower_index]
    iterp_val = vals[lower_index] + (val_step / grid_step)* (x - grid[lower_index])

    return iterp_val
end



    











In [8]:



    


grid3 = linspace(-3pi, 3pi, 19)
val3 = cos.(grid3)

cos_interp = MyInterp(grid3, val3)

@printf "cos(0) の 真値: %.5f, 補間値: %.5f\n" cos(0) cos_interp(0)
@printf "cos(π/4) の 真値: %.5f, 補間値: %.5f\n" cos(pi/4) cos_interp(pi/4)
@printf "cos(π/3) の 真値: %.5f, 補間値: %.5f\n" cos(pi/3) cos_interp(pi/3)

println("\n参考: 1/√2 = ", 1/1.4142)



    


















    






cos(0) の 真値: 1.00000, 補間値: 1.00000
cos(π/4) の 真値: 0.70711, 補間値: 0.62500
cos(π/3) の 真値: 0.50000, 補間値: 0.50000

参考: 1/√2 = 0.7071135624381276

Vectorへの対応

gridの範囲外の場合は, Noneを返すようにする





In [9]:



    


immutable MyInterp2
    grid::Array
    vals::Array
end

function (points::MyInterp2)(x::Real)
    grid = points.grid; vals = points.vals
    lower_index = searchsortedlast(points.grid, x)
    upper_index = lower_index + 1
    if lower_index == 0 || lower_index == length(grid)
        return NaN
    end

    grid_step = grid[upper_index] - grid[lower_index]
    val_step = vals[upper_index] - vals[lower_index]
    iterp_val = vals[lower_index] + (val_step / grid_step)* (x - grid[lower_index])

    return iterp_val
end



    











In [10]:



    


grid3 = linspace(-3pi, 3pi, 19)
val3 = cos.(grid3)

cos_interp = MyInterp2(grid3, val3)

@printf "cos(0) の 真値: %.5f, 補間値: %.5f\n" cos(0) cos_interp(0)
@printf "cos(π/4) の 真値: %.5f, 補間値: %.5f\n" cos(pi/4) cos_interp(pi/4)
@printf "cos(π/3) の 真値: %.5f, 補間値: %.5f\n" cos(pi/3) cos_interp(pi/3)
@printf "cos(-4π) の 真値: %.5f, 補間値: %.5f\n" cos(-4pi) cos_interp(-4pi)

println( "\nVector version: ", cos_interp.([-3pi, 2pi, 4pi, pi/4, -4pi, pi/3, 0, 2.4]) )



    


















    






cos(0) の 真値: 1.00000, 補間値: 1.00000
cos(π/4) の 真値: 0.70711, 補間値: 0.62500
cos(π/3) の 真値: 0.50000, 補間値: 0.50000
cos(-4π) の 真値: 1.00000, 補間値: NaN

Vector version: [-1.0, 1.0, NaN, 0.625, NaN, 0.5, 1.0, -0.645916]

















In [ ]:

Content source: myuuuuun/oyama_seminar2016

Similar notebooks:

notebook.community | gallery | about