14.10.2018 20:35
Yermack 15 2300 Источник


Продолжаем знакомство с очень молодым, но невероятно красивым и мощным языком программирования Julia. Шестилетняя бета наконец-таки закончилась, так что теперь можно не бояться изменений синтаксиса. И пока все спорят, хорошо или плохо начинать индексацию с единицы, взбудораженное сообщество активно закопошилось: выходят новые библиотеки, старые обновляются, стартуют серьёзные проекты, и в университетах этому языку активно учат студентов. Так не будем же отставать! Завариваем чай покрепче, потому что этой ночью будем кодить!


Подготовка к работе


Здесь есть небольшой обзор на русском, так же на хабре имеется знакомство с языком и руководство по установке. Опять же, заостряю внимание на необходимости наличия Windows Management Framework, не то будут проблемы с загрузкой пакетов.


В комплектацию JuliaPRO после обновления теперь входит только Juno. Но лично мне больше нравится Jupyter: с ним не было проблем на ноутбуке, плюс удобно работать в браузере и тут же создавать заметки и формулы, в общем, идеально для создания отчетов, слайдов или методичек.


Женим Jupyter и Julia 1.0.1
  • Скачиваем последнюю версию Юлии с официального сайта
  • Ссылка на юпитер в комплекте Анаконды дана выше, я же использовал тот, что был в старой JuliaPRO
  • Запускаем Julia. Уже можно полноценно пользоваться языком, но только в режиме интерпретатора. Выполняем команды: 
    julia>]
pkg> add IJulia
pkg> build IJulia # если не произошел автоматический build
  • Теперь в Jupyter доступно создание файла Julia 1.0.1

Для Julia существует несколько пакетов, самые же успешные из них входят в Plots в виде бэкэндов. Plots.jl — метаязык построения графика: то есть интерфейс для различных библиотек графиков. Таким образом Plots.jl на самом деле просто интерпретирует ваши команды, а затем создает графики с использованием какой-либо библиотеки графиков. Эти фоновые графические библиотеки называются бэкэндами. Самое приятное состоит в том, что вы можете использовать множество разных графических библиотек с синтаксисом Plots.jl, и мы также увидим, что Plots.jl добавляет новые функции в каждую из этих библиотек!


Устанавливаем графические пакеты

Для установки пакетов выполните команды в REPL, Juno или Jupyter:


# Pkg.add("Plots") # так добавлили пакеты до версии 0.7.0
julia>]
pkg>add Plots
pkg>add GR
pkg>add PyPlot
pkg>add Gadfly
pkg>add PlotlyJS
pkg>add UnicodePlots

Не обязательно устанавливать все пакеты, но стоит знать, что у каждого из них есть свои особенности. Я предпочитаю plotlyjs(): хоть он и не отличается быстродействием, зато очень интерактивный. Есть зум, перемещение по плоскости, а также возможность сохранения файла, причем если сохранить документ Jupyter как html, все возможности сохранятся. Так что можно добавить на сайт или сделать интерактивную презентацию. Больше информации на страницах: Plots, Gadfly


Бесконечный узор на основе простых чисел


Реализована идея статьи на хабре. В нескольких словах: что если брать координату точки и между абсциссой и ординатой применять какую-нибудь операцию, скажем, XOR или побитовое AND, а затем проверять число на простоту или на принадлежность к числам Фибоначчи, и при положительном ответе закрашивать точку в один цвет, а при отрицательном в другой? Проверим:


Для операции %
using Plots
plotlyjs()

function eratosphen(n, lst) # решето Эратосфена
    ar = [i for i=1:n]
    ar[1] = 0
    for i = 1:n
        if ar[i] != 0
            push!(lst, ar[i])
            for j = i:i:n
                ar[j] = 0
            end
        end
    end
end

ertsfn = []
eratosphen(1000, ertsfn)
# print(ertsfn)
# print( size(ertsfn) ) # -> 168
N = 80
M = 80

W1 = [in( x % y, ertsfn) for x = 1:N, y = 1:M];
W2 = [x % y for x = 1:N, y = 1:M];
p1 = spy(W1, title = "x % y is prime?")
p2 = spy(W2, title = "x % y")
plot(p1, p2, layout=(2),legend=false)


Для операции +
W1 = [in( x + y, ertsfn) for x = 1:N, y = 1:M];
W2 = [x + y for x = 1:N, y = 1:M];
p1 = spy(W1, title = "x + y is prime?")
p2 = spy(W2, title = "x + y")
plot(p1, p2, layout=(2),legend=false)


Для операции |
W1 = [in( x | y, ertsfn) for x = 1:N, y = 1:M];
W2 = [x | y for x = 1:N, y = 1:M];
p1 = spy(W1, title = "x | y is prime?")
p2 = spy(W2, title = "x | y")
plot(p1, p2, layout=(2),legend=false)


Для операции &
W1 = [in( x & y, ertsfn) for x = 1:N, y = 1:M];
W2 = [x & y for x = 1:N, y = 1:M];
p1 = spy(W1, title = "x & y is prime?")
p2 = spy(W2, title = "x & y")
plot(p1, p2, layout=(2),legend=false)


Для операции xor
W1 = [in( xor(x, y), ertsfn) for x = 1:N, y = 1:M];
W2 = [xor(x, y) for x = 1:N, y = 1:M];
p1 = spy(W1, title = "x xor y is prime?")
p2 = spy(W2, title = "x xor y")
plot(p1, p2, layout=(2),legend=false)


А теперь всё то же, но для чисел Фибонначи

Можно как обычно составлять ряд чем-то вроде:


function fib(n)
    a = 0
    b = 1
    for i = 1:n
        a, b = b, a + b
    end
    return a
end
fbncc = fib.( [i for i=1:10] )

Но воспользуемся-ка матричным представлением (подробнее):


matr_fib = n -> [1 1; 1 0]^(n-1) #  анонимная функция, возводящая матрицу в степень n-1
mfbnc = [ matr_fib( i )[1,1] for i=1:17]; # элемент 1,1 этой матрицы есть n-й элемент множества Фибоначчи
N = 100
M = N

W1 = [in( x % y, mfbnc) for x = 1:N, y = 1:M];
W2 = [in( x | y, mfbnc) for x = 1:N, y = 1:M];
p1 = spy(W1, title = "x % y ? fibonacci?")
p2 = spy(W2, title = "x | y ? fibonacci?")
plot(p1, p2, layout=(2),legend=false)


W1 = [in( xor(x, y), mfbnc) for x = 1:N, y = 1:M];
W2 = [in( x & y, mfbnc) for x = 1:N, y = 1:M];
p1 = spy(W1, title = "x xor y ? fibonacci?")
p2 = spy(W2, title = "x & y ? fibonacci?")
plot(p1, p2, layout=(2),legend=false)


Отражение фигуры относительно прямой


Данный вид отображения определяется матрицей:



где [T'], [R] и [R'] соответственно матрицы перемещения, поворота и отражения. Как это работает? Разберем на примере смещения — для точки с координатами (х, у) смещение на m по иксу и на n по игреку буде определяться преобразованием:



Эти матрицы позволяют проводить преобразования для различных многоугольников, главное записывать друг под дружкой координаты и не забывать про столбец единиц в конце. Таким образом [T]:


  • смещает прямую в началу координат вместе с преобразуемым многоугольником
  • поворачивает её до совпадения с осью Х
  • отражает все точки многоугольника относительно Х
  • затем следует обратный поворот и перенос

Более детально тема раскрыта в книге Роджерс Д., Адамс Дж. Математические основы машинной графики


А теперь закодим это!


using Plots
plotlyjs()
f = x -> 0.4x + 2 # функция задающая прямую
# иксы и игреки транспонируем в столбцы
X = [2 4 2 2]'
Y = [4 6 6 4]'
xs = [-2; 7] # отрезок принадлежащий прямой
ys = f(xs)
inptmtrx = [ X Y ones( size(X, 1), 1 ) ] # приписываем столбец единиц

m = 0
n = -f(0) # смещение по Y
displacement = [1 0 0; 0 1 0; m n 1]

a = (ys[2]-ys[1]) / (xs[2]-xs[1]) # тангенс угла наклона прямой
? = -atan(a)
rotation = [cos(?) sin(?) 0; -sin(?) cos(?) 0; 0 0 1]
reflection = [1 0 0; 0 -1 0; 0 0 1]

T = displacement * rotation * reflection * rotation^(-1) * displacement^(-1) # полная матрица преобразования
outptmtrx = inptmtrx * T
plot( X, Y)
plot!( xs, ys )
plot!( outptmtrx[:,1], outptmtrx[:,2] )


Занимательный факт: если избавиться от греческих символов и заменить первую строку на


function y=f(x,t)
    y=0.4*x + 2
endfunction,

скобки [ ] обрамляющие индексы массивов на ( ), а плоты на plot( X, Y, xs, ys, trianglenew(:,1), trianglenew(:,2) ), то данный код вполне-таки запускается в Scilab.


Работа с трехмерной графикой


Некоторые из перечисленных пакетов поддерживают построение трехмерных графиков. Но отдельно хотелось бы отметить довольно мощное средство визуализации Makie, работающий в купе с пакетами GLFW и GLAbstraction реализующими возможности OpenGL в julia. Больше информации о Makie. Его апробацию спрячем под


спойлер
using Makie

N = 51
x = linspace(-2, 2, N)
y = x
z = (-x .* exp.(-x .^ 2 .- (y') .^ 2)) .* 4

scene = wireframe(x, y, z)
xm, ym, zm = minimum(scene.limits[])
scene = surface!(scene, x, y, z)
contour!(scene, x, y, z, levels = 15, linewidth = 2, transformation = (:xy, zm))
scene


wireframe(Makie.loadasset("cat.obj"))


using FileIO

scene = Scene(resolution = (500, 500))
catmesh = FileIO.load(Makie.assetpath("cat.obj"), GLNormalUVMesh)
mesh(catmesh, color = Makie.loadasset("diffusemap.tga"))


x = Makie.loadasset("cat.obj")
mesh(x, color = :black)
pos = map(x.vertices, x.normals) do p, n
    p => p .+ (normalize(n) .* 0.05f0)
end
linesegments!(pos, color = :blue)


На этом с графикой всё. Интерактивность, анимация, 3d, большие данные или быстрое построение простейших графиков — постоянно развивающиеся пакеты удовлетворят почти любой вкус и потребности, к тому же, всё весьма легко осваивается. Смело качаем практические задания и продолжает постигать Julia!


UPD: Все вышепредставленные листинги выполнены в Jupyter с julia 0.6.4. К сожалению некоторые функции метапакета Plots либо убрали либо переименовали, так что продолжаем следить за обновлениями, а пока на скорую руку spy вполне заменяется:


julia> using GR
julia> Z = [x | y for x = 1:40, y = 1:40];
julia> heatmap(Z)

Комментарии (15)


  1. Yermack Автор
    15.10.2018 00:02
    #19235307

    outptmtrx = triangle * T
    следует заменить на 
     outptmtrx = inptmtrx * T


  1. iroln
    15.10.2018 00:22
    #19235335
    +1

    Шестилетняя бета наконец-таки закончилась, так что теперь можно не бояться изменений синтаксиса. И пока все спорят, хорошо или плохо начинать индексацию с единицы, взбудораженное сообщество активно закопошилось: выходят новые библиотеки, старые обновляются, стартуют серьёзные проекты.

    Можно пару примеров серьёзных проектов на Julia? Я наблюдал за этим языком года 4, и мне всё время казалось, что это чисто академический проект и никто не использует его в production. Если сравнить популярность Julia и Rust, то Rust за эти годы достиг гораздо больше популярности и признания и рассматривается уже как серьёзный конкурент C++. А Julia как будто не нашла свою нишу для реального применения, хотя позиционируется как быстрый и эффективный язык для параллельных вычислений, научного и инженерного применения. Если я не прав и у Джулии всё хорошо, буду только рад за них. :)


    1. eugenk
      15.10.2018 05:43
      #19235589

      Если сравнить популярность Julia и Rust, то Rust за эти годы достиг гораздо больше популярности и признания и рассматривается уже как серьёзный конкурент C++.

      Как бы Вы сравниваете несравнимое. Julia всё-таки сильно нишевой язык. Причём в достаточно консервативной области. Он просто не может расти теми же темпами что Rust! Мне кажется уже тот факт, что Julia признана на kaggle говорит о её успехе! Насчет серьёзных проектов не знаю. Но количество доступных пакетов уже радует. Причем растёт. Вобщем на мой взгляд язык очень перспективный.


      1. iroln
        15.10.2018 11:39
        #19236635

        А как на счет инструментария? IDE, интерактивный отладчик, рефакторинги, автодополнение, language server? Видел на гитхабе разные пакеты, но официальной поддержки со стороны команды разработчиков языка всего этого нет? Juno IDE выглядит очень сыро (на полноценную IDE не тянет) и он, как я понял, на основе Atom, значит будет тормозить всегда и везде.


        1. eugenk
          15.10.2018 12:19
          #19236803
          +1

          Есть на проект IDE на eclipse. Правда не знаю живой ли, интересовался этим давно. Есть на vs code. Но практика показывает, что наличие совсем уж крутых IDE народ обычно не останавливает. Я вон на яве в своё время начинал писать черт знает на чём. И ничего, потом всё появилось. Сейчас наверное вообще самый обеспеченный инструментарием язык. Так что не всё так плохо. Главное выкатили более не менее стабильную версию. Сейчас подтянут к ней пакеты, и тогда уже можно будет смотреть на дальнейшее.


      1. jamm1985
        16.10.2018 10:20
        #19240405

        А ещё в Julia можно использовать библиотеки С и Fortran. И эта возможность доступна, как говорят, «из коробки» (https://docs.julialang.org/en/v1/manual/calling-c-and-fortran-code/).


    1. Yermack Автор
      15.10.2018 11:11
      #19236491

      На оф-сайте juliacomputing есть про серьезные проекты


    1. jamm1985
      16.10.2018 10:24
      #19240433
      +1

      И ещё к преимуществам языка можно отнести наличие параллельных вычислений docs.julialang.org/en/v1/manual/parallel-computing. Можно всё через SSH делать, это очень удобно.


  1. eugenk
    15.10.2018 10:26
    #19236227

    Народ, не знаю, может кто меня поправит, но такое впечатление, что версия 1.01 (которую я сейчас поставил), стала значительно сильнее тормозить при компиляции. До того последнее с чем игрался была 0.5. Сейчас разбираться с этим особо некогда (рабочие языки в проекте verilog, java и С++), так, сел немного отвлечься и поиграться. Но в будущем планирую плотненько на джулию присесть.
    И ещё, вопрос к автору. Пробовал Ваши примеры с графикой. На spy(W1, title = «x xor y is prime?») ругается. Говорит нет такого метода. Опять же в доках не разбирался, но Вы в аргументах ничего не напутали?


    1. Yermack Автор
      15.10.2018 11:23
      #19236545

      Если про Juno, то да он довольно тормознутый. Еще довольно длительны закачка или подключение больших пакетов, поэтому, если хочется скорости, то лучше выполнять код в REPL, например загрузку пакетов. И на счет отсутствующих методов — как сказано в начале статьи, нужно скачать или обновить пакет Plots или какой-нибудь другой из перечисленных пакетов, spy у них всех работает одинаково. В качестве альтернативы можно еще попробовать Gaston

       ]add Gaston
using Gaston
...
imagesc(W1) # эквивалентно spy


      1. eugenk
        15.10.2018 12:08
        #19236773

        Да нет, я всё делаю в jupyter. Нравится мне он. Тормозит именно при компиляции. Пробовал код перезапускать — всё летает, как и должно быть. Впрочем возможен вариант, что у меня система просто уходила в своп. Довольно много всего открыто, а джулия отожрала полтора гигабайта оперативки. Надо будет попытаться поиграться на свежезагруженной.

        Насчёт spy. Всё проверил ещё раз. Не работает. В пакете Plots оно присутствует. Но очевидно не нравятся аргументы. Самое паршивое, по "?" не находится документация. Ругается вот так:
        p1 = spy(W1, title = «x % y is prime?»)
        — MethodError: no method matching findnz(::Array{Bool,2})
        Closest candidates are:
        findnz(!Matched::SparseArrays.SparseMatrixCSC{Tv,Ti}) where {Tv, Ti} at /buildworker/worker/package_linux64/build/usr/share/julia/stdlib/v1.0/SparseArrays/src/sparsematrix.jl:1337
        findnz(!Matched::SparseArrays.SparseVector{Tv,Ti}) where {Tv, Ti} at /buildworker/worker/package_linux64/build/usr/share/julia/stdlib/v1.0/SparseArrays/src/sparsevector.jl:723

        Stacktrace:
        [1] macro expansion at /home/eugenk/.julia/packages/Plots/7o1Vu/src/recipes.jl:1034 [inlined]
        [2] apply_recipe(::Dict{Symbol,Any}, ::Type{Val{:spy}}, ::UnitRange{Int64}, ::UnitRange{Int64}, ::Surface{Array{Bool,2}}) at /home/eugenk/.julia/packages/RecipesBase/Uz5AO/src/RecipesBase.jl:275
        [3] _process_seriesrecipe(::Plots.Plot{Plots.PlotlyJSBackend}, ::Dict{Symbol,Any}) at /home/eugenk/.julia/packages/Plots/7o1Vu/src/pipeline.jl:408
        [4] macro expansion at ./logging.jl:320 [inlined]
        [5] _plot!(::Plots.Plot{Plots.PlotlyJSBackend}, ::Dict{Symbol,Any}, ::Tuple{Plots.Spy}) at /home/eugenk/.julia/packages/Plots/7o1Vu/src/plot.jl:171
        [6] #plot#132(::Base.Iterators.Pairs{Symbol,String,Tuple{Symbol},NamedTuple{(:title,),Tuple{String}}}, ::Function, ::Plots.Spy) at /home/eugenk/.julia/packages/Plots/7o1Vu/src/plot.jl:57
        [7] #plot at ./none:0 [inlined]
        [8] #spy#208 at /home/eugenk/.julia/packages/RecipesBase/Uz5AO/src/RecipesBase.jl:350 [inlined]
        [9] (::getfield(Plots, Symbol("#kw##spy")))(::NamedTuple{(:title,),Tuple{String}}, ::typeof(spy), ::Array{Bool,2}) at ./none:0
        [10] top-level scope at In[19]:1


  1. eugenk
    15.10.2018 13:01
    #19236975
    +1

    P.S. Автору. Ещё один баг. Нет функции linspace(Ваш пример с 3D-графиком на Makie). Она есть в документации по 0.6.2, но в доках 1.0.1 её нет. Похоже уже убрали. Наверно что-то с версиями напутано. Потому и spy у меня не работает. Вы откуда джулию брали? Я несколько часов назад ставил с julialang.org/downloads (Generic Linux Binaries for x86 64-bit). У Вас похоже что-то более старое. В принципе этого и следовало ожидать. Недавно обновились. Пакеты ещё не стабильны. Это жизнь :))


    1. Yermack Автор
      15.10.2018 13:18
      #19237027

      У меня был юпитер от juliapro 0.6.4, все проги набил еще тогда. А потом я старый юпитер поженил с джулией 1.0.1, перекачав все пакеты. Сейчас запустил в Repl (julia 1.0.1) — spy и всё остальное работает, значит проблема в дружбе юпитера с графическими пакетами. Нда уж, с пакетами действительно придется ждать, пока их до конца перенесут на новую версию. Да и юпитер в анаконде весит много, так что меня расстроило его исключение из Джулия про. В общем со свежатиной всегда нужно быть готовым напарываться на такие ошибки, лишь бы энтузиазм не отбило


      1. eugenk
        15.10.2018 13:32
        #19237065

        У меня был юпитер от juliapro 0.6.4, все проги набил еще тогда. А потом я старый юпитер поженил с джулией 1.0.1, перекачав все пакеты.

        Тогда как говорится суду всё ясно, можно расстреливать :)))
        У Вас скорее всего в юпитере запускается всё та же 0.6.4. Я у себя джулию снёс полностью, включая каталог ~/.julia и поставил 1.0.1 начисто. Во всяком случае если у Вас нормально вызывается linspace, то это точно не 1.0.1. Можете глянуть и убедиться docs.julialang.org/en/v1/base/arrays linspace там нет.
        Вобщем-то я к этому отношусь спокойно. Такие вещи неизбежны. Для меня гораздо важнее, что с джулией я не связан по рукам и ногам библиотеками (как с питоном и numpy к примеру).


        1. Yermack Автор
          15.10.2018 14:36
          #19237357

          Да нет, в моём юпитер вполне создаются Notebook Julia 1.0.1, а на счет Макие, примеры с их сайта, а на 1.0.1 качается но не запускается. Пока для себя решил пользовать базовые возможности, а для всяких наворотов придется ждать