Hlavní navigace

Programovací jazyk Julia: pole, vektory, matice a lineární algebra

16. 6. 2016
Doba čtení: 16 minut

Sdílet

Konečně se začneme zabývat tématem souvisejícím s orientací jazyka Julia na řešení problémů z numerické matematiky. Jedná se o podporu pro práci s poli, vektory, maticemi a funkcemi dostupnými přes knihovnu LAPACK.

Obsah

1. Programovací jazyk Julia: pole, vektory, matice a lineární algebra

2. Jednorozměrná pole a dvourozměrná pole s jedním řádkem

3. Vícerozměrná pole, sloupcové a řádkové vektory

4. Vytvoření pole konstruktorem Array

5. Další konstruktory polí

6. Přístup k prvkům polí

7. Funkce range, linspace a collect

8. Funkce reshape

9. Operace slicedim pro řezy maticemi, převod matice na vektor, přístup k prvkům na diagonále

10. Lineární algebra: základní operace

11. Odkazy na Internetu

1. Programovací jazyk Julia: pole, vektory, matice a lineární algebra

Již v předchozích článcích věnovaných programovacím jazykům APL a J jsme se zmínili o tom, že jednou poměrně rozsáhlou oblastí v informatice je zpracování vektorů a matic, protože s těmito důležitými datovými strukturami se můžeme setkat v různých disciplínách, například ve finančnictví, pojišťovnictví, statistice, zpracování rozsáhlých numerických dat, simulacích, zpracování 1D a 2D signálů atd. Současně se jedná i o velmi zajímavou oblast, neboť právě kvůli co nejrychlejší práci s velkými maticemi byly vytvořeny speciální výpočetní bloky v některých superpočítačích (příkladem mohou být superpočítače Cray). Současné knihovny dokážou v případě potřeby využít jak některé rozšíření instrukčních sad (SIMD instrukce typu SSE, původně též MMX či 3DNow!), tak i programovatelné grafické akcelerátory (GPU). Práce s vektory a maticemi byla (a samozřejmě doposud je) podporována především v překladačích programovacího jazyka FORTRAN.

Překladače FORTRANu začaly být po vzniku superpočítačů vybavovány algoritmy, které dokázaly převést některé typy programových smyček na „vektorové operace“. Paralelně vznikly i specializované jazyky určené téměř výhradně pro práci s vektory i maticemi – příkladem je již zmíněná dvojice APL a J. Dnes se seznámíme s podporou práce s vektory a maticemi v programovacím jazyce Julia. Kromě toho si řekneme, jak tento jazyk dokáže využívat optimalizované a odladěné funkce dostupné ve známé knihovně LAPACK (Linear Algebra Package). Díky použití této knihovny tak může jazyk Julia a jeho uživatelé využít jednu z nejznámějších sad algoritmů pro zpracování vektorů a především matic (samotná knihovna LAPACK je překládána Fortranem 90, původní verze pak FORTRANem 77).

2. Jednorozměrná pole a dvourozměrná pole s jedním řádkem

Základní homogenní datovou strukturou, kterou programovací jazyk Julia svým uživatelům nabízí, jsou jednorozměrná pole. Všechny prvky pole mají stejný typ (ostatně právě proto je to homogenní datová struktura) a ke každému prvku je možné přistupovat přes jeho index, přičemž indexování prvků má konstantní složitost (nezáleží tedy na délce pole). Prvky v běžných jednorozměrných polích je možné měnit, takže pole jsou měnitelné datové struktury (mutable). Podívejme se nyní na způsob vytvoření jednorozměrných polí:

julia> a=[1, 2, 3, 4, 5]
5-element Array{Int64,1}:
 1
 2
 3
 4
 5

Při konstrukci pole se automaticky může zjistit datový typ prvků (resp. typ, který všem prvkům odpovídá po případné konverzi). Povšimněte si, jak se jazyk rozhoduje, který typ použít:

julia> a=[1, 2.1, 1//3]
3-element Array{Float64,1}:
 1.0
 2.1
 0.333333
 
julia> a=[1, 2, 1//3]
3-element Array{Rational{Int64},1}:
 1//1
 2//1
 1//3
 
julia> a=[1/0, -1/0, 0/0]
3-element Array{Float64,1}:
  Inf
 -Inf
  NaN
 
julia> a=[pi, pi]
2-element Array{Irrational{:π},1}:
 π = 3.1415926535897...
 π = 3.1415926535897...

Typ je možné specifikovat explicitně:

julia> Int8[1, 2, 3, 4, 5]
5-element Array{Int8,1}:
 1
 2
 3
 4
 5
 
julia> Float16[1, 2, 3, 4, 5]
5-element Array{Float16,1}:
 1.0
 2.0
 3.0
 4.0
 5.0

Pokud vynecháte čárky, vytvoří se ve skutečnosti dvourozměrné pole s jedním řádkem:

julia> a=[1 2 3 4 5]
1x5 Array{Int64,2}:
 1  2  3  4  5
 
julia> Float16[1 2 3 4 5]
1x5 Array{Float16,2}:
 1.0  2.0  3.0  4.0  5.0
 
julia> a=[1 2 3 4]
1x4 Array{Int64,2}:
 1  2  3  4

3. Vícerozměrná pole, sloupcové a řádkové vektory

Jak se však vytváří dvourozměrná pole? První pokus, který může vycházet ze zkušeností z jiných programovacích jazyků, není příliš úspěšný:

julia> [[1,2,3], [4,5,6]]
WARNING: [a,b] concatenation is deprecated; use [a;b] instead
 in depwarn at deprecated.jl:73
 in oldstyle_vcat_warning at ./abstractarray.jl:29
 in vect at abstractarray.jl:32
while loading no file, in expression starting on line 0
6-element Array{Int64,1}:
 1
 2
 3
 4
 5
 6

Problém představuje čárka vložená mezi oba vektory. Jedno z možných řešení může vypadat takto – vytvoříme vlastně pole složené ze dvou sloupců (povšimněte si chybějící čárky mezi vektory):

julia> [[1,2,3] [4,5,6]]
3x2 Array{Int64,2}:
 1  4
 2  5
 3  6

Pokud preferujete zápis po řádcích, lze použít tento alternativní způsob se středníkem. Je to sice poněkud neobvyklé, ale středník zde nahrazuje volání funkce hvcat() zmíněné níže:

julia> a=[1 2; 3 4]
2x2 Array{Int64,2}:
 1  2
 3  4

Pole se dvěma řádky a třemi sloupci:

julia> [1 2 3 ; 3 4 5]
2x3 Array{Int64,2}:
 1  2  3
 3  4  5

Kromě zápisu prvků pole do hranatých závorek lze pro konstrukci použít i funkce hcat („horizontal concatenate“), vcat („vertical concatenate“) a hvcat (kombinace obou možností se specifikací počtu sloupců):

julia> hcat(1,2,3,4)
1x4 Array{Int64,2}:
 1  2  3  4
julia> vcat(1,2,3,4)
4-element Array{Int64,1}:
 1
 2
 3
 4

U funkce hvcat() si povšimněte, že první parametr specifikuje počet sloupců a až po něm následují jednotlivé prvky:

julia> hvcat(2,1,2,3,4)
2x2 Array{Int64,2}:
 1  2
 3  4
julia> hvcat(1,1,2,3,4)
4x1 Array{Int64,2}:
 1
 2
 3
 4
julia> hvcat(4,1,2,3,4)
1x4 Array{Int64,2}:
 1  2  3  4

4. Vytvoření pole konstruktorem Array

Pro vytvoření pole s udáním typu prvků (ovšem bez inicializace jednotlivých prvků) slouží konstruktor nazvaný jednoduše Array. Při volání tohoto konstruktoru se nejprve ve složených závorkách specifikuje typ prvků a již běžně v kulatých závorkách pak rozměry pole v jednotlivých dimenzích:

help?> Array
search: Array SubArray BitArray DenseArray StridedArray mmap_array
 
  Array(dims)
 
  Array{T}(dims) constructs an uninitialized dense array with element type T.
  dims may be a tuple or a series of integer arguments. The syntax Array(T,
  dims) is also available, but deprecated.

Konstrukce pole o rozměrech 2×2 prvky typu Int8 (osmibitové celé číslo se znaménkem):

julia> a=Array{Int8}(2,2)
2x2 Array{Int8,2}:
 112  -26
  82  -34

Povšimněte si toho, že prvky skutečně nejsou inicializovány, ale obsahují obecně náhodné hodnoty:

julia> a=Array{Int8}(2,2)
2x2 Array{Int8,2}:
   80  126
 -123  -36
julia> a=Array{Int8}(2,2)
2x2 Array{Int8,2}:
 112  -121
  14   -33

Zde se nám čistě náhodou podařilo vytvořit pole s prvky s plovoucí řádovou čárkou s poloviční přesností, které mají hodnotu 0, ale samozřejmě se na to nelze spoléhat:

julia> a=Array{Float16}(10,10)
10x10 Array{Float16,2}:
 0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0
 0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0
 0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0
 0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0
 0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0
 0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0
 0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0
 0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0
 0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0
 0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0

Alternativní způsob zápisu s inicializací prvků:

julia> a = Array[[1,2], [3,4]]
2-element Array{Array{T,N},1}:
 [1,2]
 [3,4]

Při konstrukci pole se neprovádí automatické konverze, což souvisí se silným typovým systémem:

julia> a=Array{Int8}(1//2, 1//4)
ERROR: MethodError: `convert` has no method matching convert(::Type{Array{Int8,N}}, ::Rational{Int64}, ::Rational{Int64})
This may have arisen from a call to the constructor Array{Int8,N}(...),
since type constructors fall back to convert methods.
Closest candidates are:
  call{T}(::Type{T}, ::Any)
  convert{T,n}(::Type{Array{T,N}}, ::Array{T,n})
  convert{T,n,S}(::Type{Array{T,N}}, ::Array{S,n})
  ...
 in call at essentials.jl:57
 
julia> a=Array{Int8}(1.2, 3.4)
ERROR: MethodError: `convert` has no method matching convert(::Type{Array{Int8,N}}, ::Float64, ::Float64)
This may have arisen from a call to the constructor Array{Int8,N}(...),
since type constructors fall back to convert methods.
Closest candidates are:
  call{T}(::Type{T}, ::Any)
  convert{T,n}(::Type{Array{T,N}}, ::Array{T,n})
  convert{T,n,S}(::Type{Array{T,N}}, ::Array{S,n})
  ...
 in call at essentials.jl:57

5. Další konstruktory polí

Programovací jazyk Julia nabízí svým uživatelům i několik dalších funkcí určených pro konstrukci pole a současně i pro inicializaci jeho prvků. Pole o libovolných rozměrech s nulovými prvky se vytvoří funkcí zeros():

julia> b=zeros(5)
5-element Array{Float64,1}:
 0.0
 0.0
 0.0
 0.0
 0.0
julia> b=zeros(5,5)
5x5 Array{Float64,2}:
 0.0  0.0  0.0  0.0  0.0
 0.0  0.0  0.0  0.0  0.0
 0.0  0.0  0.0  0.0  0.0
 0.0  0.0  0.0  0.0  0.0
 0.0  0.0  0.0  0.0  0.0

Alternativně je možné při konstrukci pole specifikovat i typ prvků:

julia> c=zeros(Int8, 5, 5)
5x5 Array{Int8,2}:
 0  0  0  0  0
 0  0  0  0  0
 0  0  0  0  0
 0  0  0  0  0
 0  0  0  0  0

Podobná funkce, ovšem pro inicializaci všech prvků na jedničku, se jmenuje pochopitelně ones():

julia> o=ones(10,10)
10x10 Array{Float64,2}:
 1.0  1.0  1.0  1.0  1.0  1.0  1.0  1.0  1.0  1.0
 1.0  1.0  1.0  1.0  1.0  1.0  1.0  1.0  1.0  1.0
 1.0  1.0  1.0  1.0  1.0  1.0  1.0  1.0  1.0  1.0
 1.0  1.0  1.0  1.0  1.0  1.0  1.0  1.0  1.0  1.0
 1.0  1.0  1.0  1.0  1.0  1.0  1.0  1.0  1.0  1.0
 1.0  1.0  1.0  1.0  1.0  1.0  1.0  1.0  1.0  1.0
 1.0  1.0  1.0  1.0  1.0  1.0  1.0  1.0  1.0  1.0
 1.0  1.0  1.0  1.0  1.0  1.0  1.0  1.0  1.0  1.0
 1.0  1.0  1.0  1.0  1.0  1.0  1.0  1.0  1.0  1.0
 1.0  1.0  1.0  1.0  1.0  1.0  1.0  1.0  1.0  1.0

Alternativně je možné specifikovat i typ prvků:

julia> o=ones(Int8, 3, 3, 3)
3x3x3 Array{Int8,3}:
[:, :, 1] =
 1  1  1
 1  1  1
 1  1  1
 
[:, :, 2] =
 1  1  1
 1  1  1
 1  1  1
 
[:, :, 3] =
 1  1  1
 1  1  1
 1  1  1

Prvky pole lze vyplnit libovolnou hodnotou funkcí fill:

julia> fill(42, 3, 4)
3x4 Array{Int64,2}:
 42  42  42  42
 42  42  42  42
 42  42  42  42

Další užitečná funkce se jmenuje eye() a slouží k vytvoření jednotkové matice o zadaných rozměrech:

julia> eye(5,5)
5x5 Array{Float64,2}:
 1.0  0.0  0.0  0.0  0.0
 0.0  1.0  0.0  0.0  0.0
 0.0  0.0  1.0  0.0  0.0
 0.0  0.0  0.0  1.0  0.0
 0.0  0.0  0.0  0.0  1.0

Alternativně je možné i zde specifikovat typ prvků:

julia> eye(Int16, 5,5)
5x5 Array{Int16,2}:
 1  0  0  0  0
 0  1  0  0  0
 0  0  1  0  0
 0  0  0  1  0
 0  0  0  0  1

Funkce diagm() vytvoří matici, kde prvky na diagonále mají hodnoty specifikované objektem typu „range“:

julia> diagm(1:5)
5x5 Array{Int64,2}:
 1  0  0  0  0
 0  2  0  0  0
 0  0  3  0  0
 0  0  0  4  0
 0  0  0  0  5
julia> diagm(6:10)
5x5 Array{Int64,2}:
 6  0  0  0   0
 0  7  0  0   0
 0  0  8  0   0
 0  0  0  9   0
 0  0  0  0  10

V některých případech může být výhodné vytvořit pole naplněné náhodnými hodnotami. I zde nabízí jazyk Julia řešení:

julia> rand(4,5)
4x5 Array{Float64,2}:
 0.2724     0.119104   0.670161   0.835602   0.785017
 0.485103   0.207542   0.0146372  0.379119   0.329223
 0.972542   0.0750364  0.931973   0.744018   0.114805
 0.0806015  0.773603   0.162323   0.0592224  0.0124386

6. Přístup k prvkům polí

V této kapitole si ukážeme, jakým způsobem je možné přistupovat k prvkům polí. Nejprve si vytvoříme pole o velikosti 10×10 prvků, přičemž všechny prvky budou reprezentovány hodnotami s plovoucí řádovou čárkou s poloviční přesností (u rozsáhlejších polí může být použití tohoto možná poněkud neobvyklého datového typu výhodné, neboť se tím ušetří velké množství operační paměti):

julia> a=Array{Float16}(10,10)
10x10 Array{Float16,2}:
 0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0
 0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0
 0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0
 0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0
 0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0
 0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0
 0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0
 0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0
 0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0
 0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0

Pokusme se změnit hodnotu prvního prvku v poli:

julia> a[0,0]=42
ERROR: BoundsError: attempt to access 10x10 Array{Float16,2}:
  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0
  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0
  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0
  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0
  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0
  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0
  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0
  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0
  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0
  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0
  at index [0,0]
 in setindex! at array.jl:314

Vidíme, že se tato operace nepodařila, a to z toho důvodu, že se prvky indexují od jedničky a nikoli od nuly. To se sice může zdát poněkud neobvyklé, ovšem ve skutečnosti mnoho jazyků (dovolím si říci, že většina jazyků NEodvozených od céčka) zvolilo stejný přístup: Fortran, Mathematica, R, MATLAB, Lua atd. Správně tedy má příkaz vypadat takto:

julia> a[1,1]=42
42

Prvek pole se skutečně změnil:

julia> a
10x10 Array{Float16,2}:
 42.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0
  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0
  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0
  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0
  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0
  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0
  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0
  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0
  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0
  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0

Podívejme se nyní na složitější indexování, tentokrát vektoru:

julia> v=[1 2 3 4 10 -1]
1x6 Array{Int64,2}:
 1  2  3  4  10  -1

Přístup k prvnímu prvku:

julia> v[1]
1

Přístup k prvkům 2 až 4:

julia> v[2:4]
3-element Array{Int64,1}:
 2
 3
 4

Přístup ke všem prvkům:

julia> v[:]
6-element Array{Int64,1}:
  1
  2
  3
  4
 10
 -1

Pokud potřebujeme přistoupit k poslednímu prvku, není možné použít index –1 (to lze v jiných jazycích), ale používá se zde slovo end. To opět není nijak unikátní, podobně se toto slovo používá i v MATLABu:

julia> v[end]
-1

Kombinace předchozích způsobů – od čtvrtého prvku do konce vektoru:

julia> v[4:end]
3-element Array{Int64,1}:
  4
 10
 -1

Zajímavý je výběr sekvence libovolných prvků vektoru (pole), a to s využitím jiného vektoru obsahujícího indexy prvků. Povšimněte si nutnosti použití dvojitých hranatých závorek – vnější závorky představují operaci výběru prvků, vnitřní závorky vektor indexů:

julia> v[[1,5,6,2,5,5]]
6-element Array{Int64,1}:
  1
 10
 -1
  2
 10
 10

7. Funkce range, linspace a collect

Při konstrukci polí je možné využít (i když jen nepřímo) funkci pojmenovanou range(), u níž je zapotřebí si dát pozor na to, že její parametry se mohou odlišovat od stejnojmenných funkcí známých z jiných programovacích jazyků (někdy se uvádí pořadí start, stop, step atd.):

help?> range
search: range Range RangeIndex nzrange linrange histrange UnitRange StepRange
 
  range(start, [step], length)
 
  Construct a range by length, given a starting value and optional step
  (defaults to 1).

Funkce range() vytváří hodnotu typu „range“, tj. negeneruje přímo jednotlivé prvky:

julia> range(1,10)
1:10
 
julia> range(1,2,10)
1:2:19
 
julia> range(10,-1,10)
10:-1:1

Jak je tedy možné z hodnot typu „range“ vytvořit skutečné pole (či spíše vektor)? Použijeme další funkci pojmenovanou collect():

help?> collect
search: collect Collections
 
  collect(collection)
 
  Return an array of all items in a collection. For associative collections,
  returns Pair{KeyType, ValType}.
 
  collect(element_type, collection)
 
  Return an array of type Array{element_type,1} of all items in a collection.

Nenechme se zmýlit tím, že tato funkce jako svůj parametr vyžaduje kolekci, protože i objekt typu „range“ je v tomto pojetí kolekce. Pokud tedy budeme chtít vytvořit vektor hodnot od 1 do 10, použijeme tento (uznávám, že zbytečně dlouhý) zápis:

julia> collect(range(1,10))
10-element Array{Int64,1}:
  1
  2
  3
  4
  5
  6
  7
  8
  9
 10

Alternativně lze explicitně specifikovat typ prvků vytvářeného pole:

julia> collect(Float16, range(1,10))
10-element Array{Float16,1}:
  1.0
  2.0
  3.0
  4.0
  5.0
  6.0
  7.0
  8.0
  9.0
 10.0

Vzhledem k tomu, že objekt typu „range“ má i svůj literál, můžeme použít zápis:

julia> collect(3.2 : 8.6)
6-element Array{Float64,1}:
 3.2
 4.2
 5.2
 6.2
 7.2
 8.2

Popř. specifikovat jak počáteční hodnotu, tak i krok a celkový počet prvků (což je asi nejčastější způsob):

julia> collect(3.2 : 0.4 : 5)
5-element Array{Float64,1}:
 3.2
 3.6
 4.0
 4.4
 4.8

V některých případech nemusí být použití výše popsané funkce range() tím nejlepším řešením při vytváření vektoru obsahujícího sekvenci číselných hodnot. Typickým příkladem je sekvence generovaná s krokem 0.1, protože hodnotu 0.1 není možné formáty IEEE 754 single ani double přesně reprezentovat. Pokud je nutné vytvořit vektor s přesným počtem prvků, může se namísto range hodit funkce linspace(), které se předá počáteční hodnota, koncová hodnota a popř. i počet prvků vektoru. Použití funkce linspace() je tak ve skutečnosti velmi jednoduché:

julia> collect(linspace(1.0,100.0,12))
12-element Array{Float64,1}:
   1.0
  10.0
  19.0
  28.0
  37.0
  46.0
  55.0
  64.0
  73.0
  82.0
  91.0
 100.0

8. Funkce reshape

Další velmi důležitou funkcí, s níž se v praxi často setkáme, je funkce nazvaná reshape(), která dokáže změnit velikost matice a vhodným způsobem přeorganizovat prvky v původní matici. Této funkci se předávají dva parametry – prvním parametrem je vstupní pole (vektor, matice, …), druhým parametrem (popř. více parametry) pak specifikace tvaru výsledného pole. Podívejme se nejdříve na oficiální popis této funkce:

help?> reshape
search: reshape promote_shape
 
  reshape(A, dims)
 
  Create an array with the same data as the given array, but with different
  dimensions. An implementation for a particular type of array may choose
  whether the data is copied or shared.

Vytvořme si testovací vektor s dvanácti prvky (což je číslo dělitelné 2, 3, 4 i 6):

julia> a=[1 2 3 4 5 6 7 8 9 10 11 12]
1x12 Array{Int64,2}:
 1  2  3  4  5  6  7  8  9  10  11  12

Z tohoto vektoru pak snadno získáme matice o rozměrech 4×3, 3×4, 2×6, 6×2 či 1×12:

julia> reshape(a, 4, 3)
4x3 Array{Int64,2}:
 1  5   9
 2  6  10
 3  7  11
 4  8  12
julia> reshape(a, 3, 4)
3x4 Array{Int64,2}:
 1  4  7  10
 2  5  8  11
 3  6  9  12
julia> reshape(a, 2, 6)
2x6 Array{Int64,2}:
 1  3  5  7   9  11
 2  4  6  8  10  12
julia> reshape(a, 6, 2)
6x2 Array{Int64,2}:
 1   7
 2   8
 3   9
 4  10
 5  11
 6  12

Pozor ovšem na to, že počet prvků výsledné matice musí přesně odpovídat počtu prvků vstupní matice/vektoru. V tomto je Julia nekompromisní:

julia> reshape(a, 5, 2)
ERROR: DimensionMismatch("new dimensions (5,2) must be consistent with array size 12")
 in reshape at array.jl:135
 in reshape at abstractarray.jl:215

Vytvoření trojrozměrných polí je stejně snadné (a opět je nutné zachovat počet prvků):

julia> reshape(a, 2, 3, 2)
2x3x2 Array{Int64,3}:
[:, :, 1] =
 1  3  5
 2  4  6
 
[:, :, 2] =
 7   9  11
 8  10  12
julia> reshape(a, 2, 2, 3)
2x2x3 Array{Int64,3}:
[:, :, 1] =
 1  3
 2  4
 
[:, :, 2] =
 5  7
 6  8
 
[:, :, 3] =
  9  11
 10  12
julia> reshape(a, 3, 2, 2)
3x2x2 Array{Int64,3}:
[:, :, 1] =
 1  4
 2  5
 3  6
 
[:, :, 2] =
 7  10
 8  11
 9  12

Mimochodem – velikost matice lze získat funkcí size() vracející n-tici:

julia> size(reshape(a, 6,2))
(6,2)

9. Operace slicedim pro řezy maticemi, převod matice na vektor, přístup k prvkům na diagonále

Kromě funkce reshape() nabízí základní knihovna programovacího jazyka Julia i trojici dalších funkcí, které lze použít pro výběr prvků z matic předem definovaným způsobem. Jedná se o funkce pojmenované slicedim(), vec() a diag(). Podívejme se na jejich oficiální popis:

help?> slicedim
search: slicedim
 
  slicedim(A, d, i)
 
  Return all the data of A where the index for dimension d equals i.
  Equivalent to A[:,:,...,i,:,:,...] where i is in position d.
help?> vec
search: vec vecdot vecnorm Vector VecOrMat @vectorize_2arg @vectorize_1arg
 
  vec(Array) -> Vector
 
  Vectorize an array using column-major convention.
help?> diag
search: diag diagm diagind Diagonal isdiag spdiagm Bidiagonal blkdiag
 
  diag(M[, k])
 
  The kth diagonal of a matrix, as a vector. Use diagm to construct a diagonal
  matrix.

Před ukázkou možností těchto funkcí si vytvořme testovací matici o velikosti 3×4 prvky. Můžeme použít nám již známý objekt typu „range“ a funkce collect() společně s reshape():

julia> a=reshape(collect(1:12), 3, 4)
3x4 Array{Int64,2}:
 1  4  7  10
 2  5  8  11
 3  6  9  12

Získání jednotlivých řádků pole je snadné (první dimenze + index řádku):

julia> slicedim(a, 1, 1)
1x4 Array{Int64,2}:
 1  4  7  10
 
julia> slicedim(a, 1, 2)
1x4 Array{Int64,2}:
 2  5  8  11
 
julia> slicedim(a, 1, 3)
1x4 Array{Int64,2}:
 3  6  9  12

Podobně lze získat jednotlivé sloupce (druhá dimenze + index sloupce):

julia> slicedim(a, 2, 1)
3x1 Array{Int64,2}:
 1
 2
 3
 
julia> slicedim(a, 2, 2)
3x1 Array{Int64,2}:
 4
 5
 6
 
julia> slicedim(a, 2, 4)
3x1 Array{Int64,2}:
 10
 11
 12

Vidíme, že funkce slicedim() je až překvapivě flexibilní.

Funkce vec() je naproti tomu velmi jednoduchá až triviální:

julia> vec(a)
12-element Array{Int64,1}:
  1
  2
  3
  4
  5
  6
  7
  8
  9
 10
 11
 12

Podobně funkce diag(), kterou lze přečíst prvky na diagonále (tato funkce si dobře poradí i s maticí, která není čtvercová):

julia> diag(a)
3-element Array{Int64,1}:
 1
 5
 9

10. Lineární algebra: základní operace

V následující části tohoto článku si popíšeme některé funkce, které jsou nabízeny v knihovně LAPACK, ovšem již dnes se můžeme seznámit se základními operátory, které lze použít při práci s vektory a maticemi. Vytvořme si nejprve dva pomocné vektory:

julia> v1=[1,2,3]
3-element Array{Int64,1}:
 1
 2
 3
 
julia> v2=[2,2,2]
3-element Array{Int64,1}:
 2
 2
 2

Součet a rozdíl vektorů se provádí prvek po prvku:

julia> v1+v1
3-element Array{Int64,1}:
 2
 4
 6
 
julia> v1+v2
3-element Array{Int64,1}:
 3
 4
 5
 
julia> v1-v2
3-element Array{Int64,1}:
 -1
  0
  1

Vynásobení vektoru skalárem:

julia> v1*2
3-element Array{Int64,1}:
 2
 4
 6

Součin stylem prvek po prvku (první prvek s prvním prvkem atd.):

julia> v1 .* v2
3-element Array{Int64,1}:
 2
 4
 6

Podíl stylem prvek po prvku (první prvek s prvním prvkem atd.):

julia> v1 ./ v2
3-element Array{Float64,1}:
 0.5
 1.0
 1.5

Tečku lze spojit i s dalšími operátory – vždy to značí aplikaci operátoru na dvojici příslušných prvků vektoru či matice:

julia> v1 .^ v2
3-element Array{Int64,1}:
 1
 4
 9

Skalární součin vektorů:

julia> dot(v1,v2)
12

Vektorový součin vektorů:

julia> cross(v1, v2)
3-element Array{Int64,1}:
 -2
  4
 -2

Vektorový součin si lze nejlépe ukázat na jednotkových vektorech, které jsou na sebe kolmé (výsledkem je třetí vektor kolmý na oba vstupní vektory):

julia> cross([1,0,0], [0,1,0])
3-element Array{Int64,1}:
 0
 0
 1

Lineární operace odpovídající násobení skalárem:

julia> scale(v1, 10)
3-element Array{Int64,1}:
 10
 20
 30
 
julia> scale(v1, 1//10)
3-element Array{Rational{Int64},1}:
 1//10
 1//5
 3//10

Součet hodnot prvků vektorů:

julia> sum(v1)
6
 
julia> sum(v2)
6

Kombinaci tečka-operátor lze použít i pro relační operátory. Výsledkem je vektor či matice hodnot true/false:

CS24_early

julia> v1.<=v2
3-element BitArray{1}:
  true
  true
 false
julia> v1 .!= 2
3-element BitArray{1}:
  true
 false
  true

Aplikace funkce na prvky vektoru či matice nevyžaduje žádné speciální metody zápisu, minimálně pro funkce ze základní knihovny:

julia> sin(v1)
3-element Array{Float64,1}:
 0.841471
 0.909297
 0.14112

Mnohem větší množství funkcí najdeme v rozhraní knihovny LAPACK pro matice. Tyto funkce (a mnohé další) si představíme příště.

11. Odkazy na Internetu

  1. Array Programming
    https://en.wikipedia.org/wi­ki/Array_programming
  2. Discovering Array Languages
    http://archive.vector.org­.uk/art10008110
  3. no stinking loops – Kalothi
    http://www.nsl.com/
  4. Vector (obsahuje odkazy na články, knihy a blogy o programovacích jazycích APL, J a K)
    http://www.vector.org.uk/
  5. APL Interpreters
    http://www.vector.org.uk/?a­rea=interpreters
  6. APL_(programming_language
    http://en.wikipedia.org/wi­ki/APL_(programming_langu­age
  7. APL FAQ
    http://www.faqs.org/faqs/apl-faq/
  8. APL FAQ (nejnovější verze)
    http://home.earthlink.net/~swsir­lin/apl.faq.html
  9. A+
    http://www.aplusdev.org/
  10. APLX
    http://www.microapl.co.uk/
  11. FreeAPL
    http://www.pyr.fi/apl/index.htm
  12. J: a modern, high-level, general-purpose, high-performance programming language
    http://www.jsoftware.com/
  13. K, Kdb: an APL derivative for Solaris, Linux, Windows
    http://www.kx.com
  14. openAPL (GPL)
    http://sourceforge.net/pro­jects/openapl
  15. Parrot APL (GPL)
    http://www.parrotcode.org/
  16. Learning J (Roger Stokes)
    http://www.jsoftware.com/hel­p/learning/contents.htm
  17. Rosetta Code
    http://rosettacode.org/wiki/Main_Page
  18. Why APL
    http://www.acm.org/sigapl/whyapl.htm
  19. Introducing Julia/Functions
    https://en.wikibooks.org/wi­ki/Introducing_Julia/Functi­ons
  20. Functions (Julia documentation)
    http://docs.julialang.org/en/release-0.4/manual/functions/
  21. Evaluate binomial coefficients
    http://rosettacode.org/wi­ki/Evaluate_binomial_coef­ficients
  22. Ackermann function
    http://rosettacode.org/wi­ki/Ackermann_function
  23. Julia (front page)
    http://julialang.org/
  24. Julia – dokumentace
    http://docs.julialang.org/en/release-0.4/
  25. Julia – repositář na GitHubu
    https://github.com/JuliaLang/julia
  26. Julia (programming language)
    https://en.wikipedia.org/wi­ki/Julia_%28programming_lan­guage%29
  27. IJulia
    https://github.com/JuliaLan­g/IJulia.jl
  28. Introducing Julia
    https://en.wikibooks.org/wi­ki/Introducing_Julia
  29. Julia: the REPL
    https://en.wikibooks.org/wi­ki/Introducing_Julia/The_REPL
  30. Introducing Julia/Metaprogramming
    https://en.wikibooks.org/wi­ki/Introducing_Julia/Meta­programming
  31. Month of Julia
    https://github.com/DataWo­okie/MonthOfJulia
  32. Learn X in Y minutes (where X=Julia)
    https://learnxinyminutes.com/doc­s/julia/
  33. New Julia language seeks to be the C for scientists
    http://www.infoworld.com/ar­ticle/2616709/application-development/new-julia-language-seeks-to-be-the-c-for-scientists.html
  34. Julia: A Fast Dynamic Language for Technical Computing
    http://karpinski.org/publi­cations/2012/julia-a-fast-dynamic-language
  35. The LLVM Compiler Infrastructure
    http://llvm.org/
  36. Julia: benchmarks
    http://julialang.org/benchmarks/
  37. Type system
    https://en.wikipedia.org/wi­ki/Type_system
  38. Half-precision floating-point format
    https://en.wikipedia.org/wiki/Half-precision_floating-point_format

Byl pro vás článek přínosný?

Autor článku

Vystudoval VUT FIT a v současné době pracuje na projektech vytvářených v jazycích Python a Go.