# introduction-to-julia

## Introduction to Julia (slides show)¶

Let's start with a quick overview of the basic syntax, emphasising differences with Python.

Why Jula (slides show)

• David P. Sanders
• Department of Physics, Faculty of Sciences
• National University of Mexico (UNAM)
• GitHub: dpsanders

## Variables and arithmetic¶

Numeric values in Julia work similarly to Python:

Variables in Julia are created as in Python, with a simple assignment operation; variable names can be arbitrary unicode characters. Many may be entered in LaTeX notation, using tab substitution: type \alpha<TAB>. There is also tab completion on partial names: \alp<TAB>

In [1]:
x = 3
Out[1]:
3
In [2]:
y = 5
Out[2]:
5
In [3]:
α = 3;  = 10
Out[3]:
10

Functions use parentheses (round brackets, ()) around the arguments being passed. println prints its arguments, followed by a new line. [print omits the new line.]

In [4]:
println("α = ", α)
α = 3

Simple functions may be defined with a nice mathematical syntax; * is not needed in simple expressions:

In [5]:
f(x) = 2x^2 + 3x + 1
g(x) = f(x) - (2x+1)*(x+1)
Out[5]:
g (generic function with 1 method)
In [6]:
f(3)
Out[6]:
28
In [7]:
g(3.5)
Out[7]:
0.0

### Variable substitution¶

The values of variables may be substituted into strings in a simple way using the $operator: In [8]: # Variable substitution with$:
name = "David"
greeting = "Hello, $name" println(greeting) Hello, David More complicated expressions are wrapped in parentheses: In [9]: μ = 3 println("The sine of$μ is $(sin(μ))") The sine of 3 is 0.1411200080598672 ### Numerical types¶ There are numerical types with different precisions: typing Float<TAB> or Int<TAB> will provide a list. Currently, in arithmetic calculations types are promoted to the machine type. (This looks likely to change soon.) Machine integers! In [10]: a = int(1e16) a * 10 Out[10]: 100000000000000000 In [11]: a = int8(1) b = int8(2) a + b Out[11]: 3 In [12]: typeof(ans) # ans is the last result Out[12]: Int8 These promotion rules are defined in int.jl. ### Arbitrary-precision arithmetic¶ Arbitrary-precision integers and floating points are available through the types BigInt and BigFloat. The function big converts a number into the corresponding Big type: In [13]: big(10) Out[13]: 10 In [14]: typeof(ans) Out[14]: BigInt (constructor with 10 methods) Note that, unlike in Python, integers are not automatically promoted to arbitrary-precision integers. Exercise: Calculate powers of 10 using standard integers and BigInts In [15]: 10^18 Out[15]: 1000000000000000000 In [16]: 10**5 syntax: use "^" instead of "**" while loading In[16], in expression starting on line 1 In [17]: 10^19 Out[17]: -8446744073709551616 In [18]: ten = big(10) Out[18]: 10 In [19]: ten^109 Out[19]: 10000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000 In [20]: i = int8(10) Out[20]: 10 In [21]: i * big(11) Out[21]: 110 In [22]: typeof(ans) Out[22]: BigInt (constructor with 10 methods) In [23]: typemax(Int64) Out[23]: 9223372036854775807 In [24]: 2^63-1 Out[24]: 9223372036854775807 In [25]: typemin(Int64) Out[25]: -9223372036854775808 In [26]: -2^63 Out[26]: -9223372036854775808 ### Complex numbers¶ Complex numbers are written using im for the imaginary part: In [27]: a = 7 c = (1+3.5im) * a Out[27]: 7.0 + 24.5im In [28]: c.im Out[28]: 24.5 In [29]: c.re, c.im Out[29]: (7.0,24.5) In [30]: c * conj(c) # conj is a function that returns the conjugate of a complex number Out[30]: 649.25 + 0.0im (Tuples behave similarly to Python.) ### Rational numbers¶ Rational numbers are also built into Julia; they are created using the // operator: In [31]: 3//4 Out[31]: 3//4 In [32]: typeof(ans) Out[32]: Rational{Int64} (constructor with 1 method) In [33]: (big(3)//4)^50 Out[33]: 717897987691852588770249//1267650600228229401496703205376 In [34]: typeof(ans) Out[34]: Rational{BigInt} (constructor with 1 method) In [35]: 3//4 + 5//6 Out[35]: 19//12 Operators are a convenient way of writing functions: In [36]: +(3, 4) Out[36]: 7 In [37]: //(3, 4) Out[37]: 3//4 In [38]: // Out[38]: // (generic function with 8 methods) We see that // is a function, implemented as a series of methods. We can see what these methods are: In [39]: methods(//) Out[39]: 8 methods for generic function //: In [40]: 3 + -34 Out[40]: -31 The expression n::Integer is a type annotation that specifies that the method applies when its first argument is of type Integer. Clicking on the file name takes us directly to the Julia standard library source code on GitHub where these functions are defined! ## Vectors: an "equivalent" of Python lists andnumpy arrays¶ To store several values "in one variable", we can try to imitate using a "list" as we would in Python: In [41]: l = [3, 4, 5] Out[41]: 3-element Array{Int64,1}: 3 4 5 In [42]: typeof(l) Out[42]: Array{Int64,1} In Julia these objects are called Arrays. The curly braces indicate type parameters of the Array type. The first is the type of element contained in the Array (all must be of the same type) and the second the number of dimensions. Exercise: Try to create an array in this way with elements of different types. What happens? Exercise: What does the following syntax do?: l = {3, 4, 7.5} In [43]: l = [3, 4, 7.5] Out[43]: 3-element Array{Float64,1}: 3.0 4.0 7.5 In [44]: l = [3., 4, 5] Out[44]: 3-element Array{Float64,1}: 3.0 4.0 5.0 In [45]: l = [3., "a"] Out[45]: 2-element Array{Any,1}: 3.0 "a" In [46]: l = [3., 'a'] Out[46]: 2-element Array{Float64,1}: 3.0 97.0 In [47]: l = {3., 4, "hello", [3, 4]} Out[47]: 4-element Array{Any,1}: 3.0 4 "hello" [3,4] ### Indexing¶ The indices of Julia arrays are numbered starting at 1, unlike Python (where they are numbered starting at 0). In [48]: l[1] Out[48]: 3.0 The syntax for ranges is similar to that for Python: In [49]: l[1:2] Out[49]: 2-element Array{Any,1}: 3.0 4 However, the limits must be explicitly specified: In [50]: l[2:end] # Use end explicitly Out[50]: 3-element Array{Any,1}: 4 "hello" [3,4] In [51]: l[1:end-1] Out[51]: 3-element Array{Any,1}: 3.0 4 "hello" In [52]: l[-1] BoundsError() while loading In[52], in expression starting on line 1 in getindex_3B_2413 at /home/tj2/Documents/julia/usr/bin/../lib/julia/sys.so Julia Arrays, like Python lists, but unlike numpy arrays, are dynamic. However, the syntax is rather different from Python -- to add an element at the end of the list, we write In [53]: l = [3,4,5] l + l Out[53]: 3-element Array{Int64,1}: 6 8 10 In [54]: names(l) Out[54]: 0-element Array{Any,1} In [55]: l = [3, 4, 5] push!(l, 7) Out[55]: 4-element Array{Int64,1}: 3 4 5 7 In [56]: methods(sizehint) Out[56]: 11 methods for generic function sizehint: In [57]: push! Out[57]: push! (generic function with 22 methods) In [58]: methods(push!) Out[58]: 22 methods for generic function push!: In [60]: l = [3, 4, 5] push!(l, 12.001) InexactError() while loading In[60], in expression starting on line 2 in push! at array.jl:451 In [61]: l = [3.0, 4, 5] push!(l, 12.001) Out[61]: 4-element Array{Float64,1}: 3.0 4.0 5.0 12.001 In [62]: 12.0 == 12 Out[62]: true In [63]: push!(l, 12.1) Out[63]: 5-element Array{Float64,1}: 3.0 4.0 5.0 12.001 12.1 In [64]: l = [3, 4, 5, 12] append!(l, [10., 11., 12]) Out[64]: 7-element Array{Int64,1}: 3 4 5 12 10 11 12 push! replaces append in Python. There are no methods of objects as in Python; rather, we use functions and send the object as an argument of the function. The exclamation mark, or bang, (!) indicates that the function modifies its argument; this is a standard convention in Julia. Arrays which have been defined with a certain type cannot acquire elements of a different type: In [65]: l = [3, 4, 5] push!(l, "hello") convert has no method matching convert(::Type{Int64}, ::ASCIIString) while loading In[65], in expression starting on line 2 in push! at array.jl:451 Arrays work as mathematical vectors, with the sum of two vectors and scalar multiplication being defined: In [66]: a = [1.1, 2.2, 3.3] b = [4.4, 5.5, 6.6] Out[66]: 3-element Array{Float64,1}: 4.4 5.5 6.6 In [67]: a + b Out[67]: 3-element Array{Float64,1}: 5.5 7.7 9.9 In [68]: 3.5 * a Out[68]: 3-element Array{Float64,1}: 3.85 7.7 11.55 However, operators are, in general, not treated in an elementwise fashion (as they would be e.g. in numpy): In [69]: a * b * has no method matching *(::Array{Float64,1}, ::Array{Float64,1}) while loading In[69], in expression starting on line 1 Rather, elementwise operations use a Matlab-like syntax, with an extra . before the symbol for the operator: In [70]: a .* b Out[70]: 3-element Array{Float64,1}: 4.84 12.1 21.78 There are many useful operations on vectors predefined, without needing to explicitly import them. In [71]: dot(a,b) Out[71]: 38.72 In [72]: ?dot dot (generic function with 7 methods) In [73]: cross(a, b) Out[73]: 3-element Array{Float64,1}: -3.63 7.26 -3.63 In [74]: ?cross cross (generic function with 1 method) In [75]: norm(a) Out[75]: 4.115823125451335 In [76]: ?norm Base.norm(A[, p]) Compute the "p"-norm of a vector or the operator norm of a matrix "A", defaulting to the "p=2"-norm. For vectors, "p" can assume any numeric value (even though not all values produce a mathematically valid vector norm). In particular, "norm(A, Inf)" returns the largest value in "abs(A)", whereas "norm(A, -Inf)" returns the smallest. For matrices, valid values of "p" are "1", "2", or "Inf". (Note that for sparse matrices, "p=2" is currently not implemented.) Use "vecnorm()" to compute the Frobenius norm. Use help or ? (before the command) to obtain help: In [77]: help(dot) dot (generic function with 7 methods) In [78]: ?dot dot (generic function with 7 methods) In [79]: transpose(a) Out[79]: 1x3 Array{Float64,2}: 1.1 2.2 3.3 In [80]: a' Out[80]: 1x3 Array{Float64,2}: 1.1 2.2 3.3 In [81]: transpose(a) == a' Out[81]: true In [82]: M = [[2,1], [1,1]] # 1-d Out[82]: 4-element Array{Int64,1}: 2 1 1 1 In [83]: push!(l, 12.0) Out[83]: 4-element Array{Int64,1}: 3 4 5 12 In [84]: M = [2. 1; 1 1] # 2-d Out[84]: 2x2 Array{Float64,2}: 2.0 1.0 1.0 1.0 In [85]: M = reshape([1:8], (2,2,2)) Out[85]: 2x2x2 Array{Int64,3}: [:, :, 1] = 1 3 2 4 [:, :, 2] = 5 7 6 8 In [86]: a b Out[86]: 38.72 In [87]: dot(a,b) == a b Out[87]: true In [88]: a × b Out[88]: 3-element Array{Float64,1}: -3.63 7.26 -3.63 In [89]: cross(a,b) == a × b Out[89]: true I used \cdot In [90]: Out[90]: dot (generic function with 7 methods) [Note that in the Julia command-line REPL, typing ? puts it immediately into a special help mode. Similarly, ; puts it into shell mode, in which commands are sent straight to the shell.] ## Control flow: drop the colon (:) and add end¶ White space in Julia is not significant. Commands on one line can be separated by ;. Blocks must finish with end In [91]: i = 0 while i < 5 print("$i\t")
i += 1
end
0	1	2	3	4
In [92]:
total = 0
for i = 1:26
total += 2^i
end
println("Sum is $total") Sum is 134217726 Here, 1:26 is a range object which may be iterated over. In [93]: typeof(1:10) Out[93]: UnitRange{Int64} (constructor with 1 method) We can construct an array from this by enclosing it in square brackets: In [94]: collect(1:10) Out[94]: 10-element Array{Int64,1}: 1 2 3 4 5 6 7 8 9 10 In [95]: ?collect Base.collect(collection) Return an array of all items in a collection. For associative collections, returns (key, value) tuples. Base.collect(element_type, collection) Return an array of type "Array{element_type,1}" of all items in a collection. In [96]: [1:2:10, 17] Out[96]: 6-element Array{Int64,1}: 1 3 5 7 9 17 In [97]: help("dot") Base.LinAlg.BLAS.dot(n, X, incx, Y, incy) Dot product of two vectors consisting of "n" elements of array "X" with stride "incx" and "n" elements of array "Y" with stride "incy". Exercise: Implement the Babylonian method for calculating the square root of a positive number$y$, via the iteration $$x_{n+1} = \textstyle \frac{1}{2} (x_n + \frac{y}{x_n})$$ In [98]: x = 7 # Babylonian method r1 = big(x) for i=1:10 r2 = (r1+x//r1)//2 r1 = r2 println("$i: $r2 is$(float(r2))")
end
sqrt(x)
1: 4//1 is 4e+00
2: 23//8 is 2.875e+00
3: 977//368 is 2.654891304347826086956521739130434782608695652173913043478260869565217391304333e+00
4: 1902497//719072 is 2.645767044190289706733122691469004494682034622402207289395220506430510435672636e+00
5: 7238946623297//2736064645568 is 2.645751311111369320474788060415261270876928274529823742547961659814056687694087e+00
6: 104804696428033056657448577//39612451854313553433195392 is 2.645751311064590590502029293938515493859933284910652657274139597282698886360655e+00
7: 21968048786744329890159858491687242656637852949560577//8303141982722814210287293056980408998691454066714368 is 2.645751311064590590501615753639260425710259215401356797324429072042618334785781e+00
8: 965190334993558048952747224707737662553546023969581717336893753948305178090306054973804711697153676048897//364807656319439656893574197813762741889129729119092813240115363343939412817707718620985339231655544540672 is 2.645751311064590590501615753639260425710259183082450180368334459201068823230298e+00
9: 1863184765529953614442676969410516768705348877426018592874619682164981370927579559773987864606096009580887559742891544147556039040815814753743851100010309178655530447444083262328821461724034060383304637370397697//704217648022349512771442189018194725339511495157190929516136319260161990261007461575776992568187313708188185637682782310132281240170939460318854442124039911242297405920199037531309472261446561026954943354477568 is 2.645751311064590590501615753639260425710259183082450180368334459201068823230264e+00
10: 6942914941005816452906820447809898528401948252945007180207641686978710969653763380353810662539739407094605662843293836754945636431635508103618963852112274901648239421682834044604657380796983687411621541616911160282084512012755071820846262756649514791264502834906610162581301591431279046989968494539999498942229056113659738389773073293970820133886115413882850641756087247746883347629340780021440048884451568277591449266177//2624175186825153359249378049906502688043810272014845844292434248994264560667710669982033339453245573538963306071655131184848403463172113026723021373675261241203122087530398269881154662495920113030441863771241633914463530742637259717641425928569498576675875622441689707127692085122998906191651056070668909263006767739390431278990942834266700554990212884068483635241604906439758589382529929519521217527464936293523650721792 is 2.645751311064590590501615753639260425710259183082450180368334459201068823230264e+00
Out[98]:
2.6457513110645907
In [99]:
= norm
Out[99]:
norm (generic function with 15 methods)
In [100]:
(a)
Out[100]:
4.115823125451335
In [101]:
snowman(x) = x^3
Out[101]:
snowman (generic function with 1 method)
In [102]:
= snowman
(3)
Out[102]:
27

### Short-circuit evaluation¶

In [103]:
a = 3
a < 5 && println("Small")   # evaluate the second statement only if the first is true;  semantics of if-then

a > 10 || println("Small")  # semantics of if not-then

# equal to

if a < 5
println("small, if-then")
end

if a > 10
println("small, if not-then")
end
Small
Small
small, if-then
In [104]:
a == 3 ? println("Hello") : println("Not true")

# equal to

if a < 5
println("small, if-then")
else
println("small, if not-then")
end
Hello
small, if-then

### Array comprehensions¶

There is an equivalent of list comprehensions in Python, as follows. Note that the array construction syntax is quite flexible.

In [105]:
squares = [i^2 for i in [1:2:10, 7]]
Out[105]:
6-element Array{Any,1}:
1
9
25
49
81
49
In [106]:
sums = [i+j for i=1:5, j=1:5] # boardcasting
Out[106]:
5x5 Array{Int64,2}:
2  3  4  5   6
3  4  5  6   7
4  5  6  7   8
5  6  7  8   9
6  7  8  9  10
In [107]:
sums = [i+j+k for i=1:5, j=1:5, k=1:5]
Out[107]:
5x5x5 Array{Int64,3}:
[:, :, 1] =
3  4  5   6   7
4  5  6   7   8
5  6  7   8   9
6  7  8   9  10
7  8  9  10  11

[:, :, 2] =
4  5   6   7   8
5  6   7   8   9
6  7   8   9  10
7  8   9  10  11
8  9  10  11  12

[:, :, 3] =
5   6   7   8   9
6   7   8   9  10
7   8   9  10  11
8   9  10  11  12
9  10  11  12  13

[:, :, 4] =
6   7   8   9  10
7   8   9  10  11
8   9  10  11  12
9  10  11  12  13
10  11  12  13  14

[:, :, 5] =
7   8   9  10  11
8   9  10  11  12
9  10  11  12  13
10  11  12  13  14
11  12  13  14  15

## Matrices¶

Square brackets with commas gives a one-dimensional vector. This is printed in a way that treats it as if it were a column vector (although there is in fact no difference between a one-dimensional row vector and column vector).

In [108]:
v = [3, 4, 5]
Out[108]:
3-element Array{Int64,1}:
3
4
5

To create explicit matrices, Matlab-style notation is used. If we omit the commas, something different happens: we now obtain a two-dimensional Array, i.e. a matrix, of size $1 \times n$. [Recall that in the standard notation for matrices, an $m \times n$ matrix has $m$ rows and $n$ columns.]

In [109]:
row_vec = [3 4 5]
Out[109]:
1x3 Array{Int64,2}:
3  4  5
In [110]:
v == row_vec
Out[110]:
false

We can also use the transpose operator, '. [This is actually the conjugate-transpose operator, which also takes the complex conjugate of complex numbers. Transpose without conjugate is denoted .']

In [111]:
row_vec = [2im, 2]'
Out[111]:
1x2 Array{Complex{Int64},2}:
0-2im  2+0im
In [112]:
row_vec = [2im, 2].'
Out[112]:
1x2 Array{Complex{Int64},2}:
0+2im  2+0im

A complete matrix may be constructed using a semicolon (;) to separate rows:

In [113]:
M = [1 2im; 3 4]
Out[113]:
2x2 Array{Complex{Int64},2}:
1+0im  0+2im
3+0im  4+0im

As in numpy, it may also be created using a reshape:

In [114]:
M = reshape([1 2 3 4], (2,2))

# equal to

N = reshape([1, 2, 3, 4], (2,2))
Out[114]:
2x2 Array{Int64,2}:
1  3
2  4
In [115]:
M == N
Out[115]:
true

Here, as in Python, (2,2) denotes an (immutable) tuple:

In [116]:
t = (2, 2)
typeof(t)
Out[116]:
(Int64,Int64)

There is an important difference in the way that Python and Julia treat slices of matrices. While in Python a one-dimensional slice in either direction returns a 1-dimensional vector, in Julia there is a difference. A vertical one-dimensional slice gives a 1-dimensional vector (a "column vector"):

In [117]:
M[:,1]
Out[117]:
2-element Array{Int64,1}:
1
2

However, a horizontal one-dimensional slice produces a $1 \times n$ matrix:

In [118]:
M[1,:]
Out[118]:
1x2 Array{Int64,2}:
1  3

This is the same result that is produced using the following Matlab-like syntax:

In [119]:
[1 2]
Out[119]:
1x2 Array{Int64,2}:
1  2

### Random numbers¶

The Mersenne Twister (pseudo-)random number generator is built-in to Julia:

In [120]:
rand()
Out[120]:
0.40186768641976944
In [121]:
rand(5)
Out[121]:
5-element Array{Float64,1}:
0.344519
0.617192
0.829609
0.981215
0.477052
In [122]:
x = rand(5, 5)
Out[122]:
5x5 Array{Float64,2}:
0.634627  0.121441   0.793247  0.78431   0.965882
0.989307  0.616908   0.353485  0.710715  0.635409
0.352762  0.118025   0.185011  0.280691  0.28924
0.748569  0.195795   0.279713  0.831955  0.888495
0.765     0.0837347  0.537341  0.429161  0.147381
In [123]:
y = rand(Int32,3,3)
Out[123]:
3x3 Array{Int32,2}:
-113324951   196613133  -1363228423
-2088630349  1597727404   -358295295
-1080816523  1898363541    571418313
In [124]:
?rand
Base.rand() -> Float64

Generate a "Float64" random number uniformly in [0,1)

Base.rand(rng::AbstractRNG[, dims...])

Generate a random "Float64" number or array of the size specified
by dims, using the specified RNG object. Currently,
"MersenneTwister" is the only available Random Number Generator
(RNG), which may be seeded using srand.

Base.rand(dims or [dims...])

Generate a random "Float64" array of the size specified by dims

Base.rand(t::Type[, dims...])

Generate a random number or array of random numbes of the given
type.

Base.rand(r[, dims...])

Pick a random element or array of random elements from range "r"
(for example, "1:n" or "0:2:10").
In [125]:
methods(rand)
Out[125]:
33 methods for generic function rand:

### Matrix multiplication¶

In [126]:
v = [1, 2]
Out[126]:
2-element Array{Int64,1}:
1
2
In [127]:
v*v
* has no method matching *(::Array{Int64,1}, ::Array{Int64,1})
In [128]:
v^2
* has no method matching *(::Array{Int64,1}, ::Array{Int64,1})

in power_by_squaring at intfuncs.jl:56
in ^ at intfuncs.jl:86
In [129]:
dot(v, v)
Out[129]:
5
In [130]:
M = [2 1; 1 1]
Out[130]:
2x2 Array{Int64,2}:
2  1
1  1
In [131]:
dot(M, v)
dot has no method matching dot(::Array{Int64,2}, ::Array{Int64,1})

Matrix multiplication uses the * operator:

In [132]:
M * v
Out[132]:
2-element Array{Int64,1}:
4
3
In [133]:
@which M*v
Out[133]:
*{T,S}(A::AbstractArray{T,2},x::AbstractArray{S,1}) at linalg/matmul.jl:71

Exercise: Use the power method to calculate the largest eigenvalue $\lambda_1$ of the matrix $M = \begin{pmatrix} 2 & 1 \\ 1 & 1 \end{pmatrix}$. In this method, we start from an arbitrary non-zero vector $\mathbf{w}$, and repeatedly apply $M$ to it, thus calculating powers of the matrix $M$ applied to $\mathbf{w}$. The resulting vector converges to the eigenvector $\mathbf{v}_1$ corresponding to $\lambda_1$.

In [134]:
w = [1., 1]
M = reshape([2., 1, 1, 1], (2,2))

M, w
Out[134]:
(
2x2 Array{Float64,2}:
2.0  1.0
1.0  1.0,

[1.0,1.0])
In [135]:
eig(M)
Out[135]:
([0.381966,2.61803],
2x2 Array{Float64,2}:
0.525731  -0.850651
-0.850651  -0.525731)
In [136]:
?eig
Base.eig(A,[irange,][vl,][vu,][permute=true,][scale=true]) -> D, V

Computes eigenvalues and eigenvectors of "A". See "eigfact()"
for details on the "balance" keyword argument.

julia> eig([1.0 0.0 0.0; 0.0 3.0 0.0; 0.0 0.0 18.0])
([1.0,3.0,18.0],
3x3 Array{Float64,2}:
1.0  0.0  0.0
0.0  1.0  0.0
0.0  0.0  1.0)

"eig" is a wrapper around "eigfact()", extracting all parts of
the factorization to a tuple; where possible, using "eigfact()"
is recommended.

Base.eig(A, B) -> D, V

Computes generalized eigenvalues and vectors of "A" with respect
to "B".

"eig" is a wrapper around "eigfact()", extracting all parts of
the factorization to a tuple; where possible, using "eigfact()"
is recommended.
In [137]:
w0 = [1., 1]
w = copy(w0)

for i in 1:10
w_new = M*w
println(w_new)
end
[3.0,2.0]
[3.0,2.0]
[3.0,2.0]
[3.0,2.0]
[3.0,2.0]
[3.0,2.0]
[3.0,2.0]
[3.0,2.0]
[3.0,2.0]
[3.0,2.0]
In [138]:
for i =1:10
t = 3
end
t
Out[138]:
3

## Linear algebra¶

Julia has built-in linear algebra, not only using LAPACK, but now also generic routines that work for arbitrary element types, implemented completely in Julia.

For example, given a matrix $A$, the LU-decomposition of $A$ is equivalent to Gaussian elimination; it expresses $A$ as the product $A = LU$, with $L$ a lower-triangular and $U$ an upper-triangular matrix.

This is implemented in pure Julia for arbitrary element types. When the elements are standard floating-point numbers, it uses the corresponding fast LAPACK implementation.

In [139]:
M = rand(100, 100)
eig(M)
Out[139]:
(Complex{Float64}[50.0183+0.0im,1.84743+2.30961im,1.84743-2.30961im,-2.72534+0.880116im,-2.72534-0.880116im,-2.86948+0.0im,2.64375+0.964703im,2.64375-0.964703im,2.79553+0.0im,0.196765+2.72857im  …  -0.867929-0.457136im,0.478782+0.882849im,0.478782-0.882849im,-0.120951+0.762254im,-0.120951-0.762254im,-0.54004+0.0im,0.471338+0.0im,0.165157+0.118936im,0.165157-0.118936im,-0.288467+0.0im],
100x100 Array{Complex{Float64},2}:
-0.101511+0.0im   -0.0115356-0.039612im    …   -0.0524037+0.0im
-0.0952353+0.0im   0.00489796-0.00599155im      -0.0682674+0.0im
-0.108918+0.0im    0.0161753-0.092386im         -0.114214+0.0im
-0.111416+0.0im     0.102416+0.0351175im        -0.103224+0.0im
-0.0977978+0.0im   -0.0198729+0.0722435im      -0.00218148+0.0im
-0.0836507+0.0im     0.082544+0.0627261im   …    -0.036751+0.0im
-0.0963353+0.0im   -0.0388372+0.0032277im        0.0954342+0.0im
-0.0929018+0.0im     0.029581+0.0575031im        0.0794094+0.0im
-0.108383+0.0im   -0.0832873-0.160539im         0.0634216+0.0im
-0.0902722+0.0im   -0.0375187+0.0593948im       -0.0071044+0.0im
-0.10374+0.0im    -0.102945+0.00231717im  …   -0.0538666+0.0im
-0.106464+0.0im    -0.104489-0.00207164im      -0.0406458+0.0im
-0.0960566+0.0im    0.0159221-0.049508im          -0.16083+0.0im
⋮                                 ⋱
-0.0946263+0.0im   -0.0375218+0.0672119im          -0.1078+0.0im
-0.0908183+0.0im  -0.00540686+0.142578im           0.12745+0.0im
-0.104126+0.0im   -0.0670472+0.0553909im   …    -0.158563+0.0im
-0.0999715+0.0im    0.0382579-0.0947829im       -0.0696019+0.0im
-0.0986926+0.0im    -0.084612-0.12682im           0.262889+0.0im
-0.0949404+0.0im   -0.0400257+0.023941im          -0.11199+0.0im
-0.107486+0.0im   -0.0344824-0.0759473im         0.087478+0.0im
-0.0944187+0.0im    -0.021001+0.0461968im   …    0.0110733+0.0im
-0.0999898+0.0im    0.0392409-0.0264316im       0.00222877+0.0im
-0.0991767+0.0im  -0.00947725+0.0518241im        0.0205667+0.0im
-0.108745+0.0im    0.0184766-0.0401366im         0.232166+0.0im
-0.0993031+0.0im    0.0796098-0.0182621im        0.0712234+0.0im)
In [140]:
M = rand(100, 100)
M2 = map(big, M)
Out[140]:
100x100 Array{BigFloat,2}:
2.225506282223690401878002376179210841655731201171875e-01  …  6.56851621473254265737296009319834411144256591796875e-01
7.319635804294619685350653526256792247295379638671875e-01     9.672824011227445506477806702605448663234710693359375e-01
5.340947103998476830355457423138432204723358154296875e-01     9.836664954220417644847884730552323162555694580078125e-01
3.162661116152143403468244287068955600261688232421875e-01     6.32452881486862406035243111546151340007781982421875e-01
2.223050262746237848432429018430411815643310546875e-01        9.4728250126541180264894137508235871791839599609375e-02
9.375878604020722040246482720249332487583160400390625e-01  …  5.216604894602110409351780617726035416126251220703125e-01
2.82637383973149969307314677280373871326446533203125e-02      7.1043767348334707634194273850880563259124755859375e-01
9.31414288307105042719058474176563322544097900390625e-01      6.843501947718078071147829177789390087127685546875e-01
9.08755575606785992448521938058547675609588623046875e-01      7.7541670796339889903947550919838249683380126953125e-02
9.90974279386318102069708402268588542938232421875e-01         3.446265247750661675496530733653344213962554931640625e-01
3.62011301348933667298979344195686280727386474609375e-01   …  5.674886600100188527306954711093567311763763427734375e-01
9.70427871977250333657138980925083160400390625e-01            1.845040577422647487537687993608415126800537109375e-01
9.0041572068816577711913851089775562286376953125e-02          2.86408155710278666816748227574862539768218994140625e-01
⋮                                                          ⋱
5.51151260187439984150614691316150128841400146484375e-01      3.306386894959143507577437048894353210926055908203125e-01
3.497623439347330442927841431810520589351654052734375e-01     8.4233377002282594503412838093936443328857421875e-01
2.454565323093833573153688121237792074680328369140625e-01  …  3.436356112336635959536579321138560771942138671875e-01
9.3207496351350638263966175145469605922698974609375e-01       7.007260372885506516382747577154077589511871337890625e-01
8.144309263423663569625432501197792589664459228515625e-01     3.333819991414610495183978855493478477001190185546875e-01
7.57929320416632634760389919392764568328857421875e-02         6.67666265761114541277265743701718747615814208984375e-01
6.366199276186090827422958682291209697723388671875e-01        3.03926419592160268479119622497819364070892333984375e-01
5.72819224497360668379997150623239576816558837890625e-01   …  6.6274981934506360659042911720462143421173095703125e-01
9.56949160642551976962977278162725269794464111328125e-01      8.29108303980233163343882551998831331729888916015625e-01
7.970001203346741558419807915925048291683197021484375e-01     1.888462575949938493380386717035435140132904052734375e-01
9.0568890895109799288320573396049439907073974609375e-01       4.966056054863334789928330792463384568691253662109375e-01
2.048115155536647247203063670895062386989593505859375e-01     3.74601565881409026559367703157477080821990966796875e-01
In [141]:
?map
Base.map(f, c...) -> collection

Transform collection "c" by applying "f" to each element. For
multiple collection arguments, apply "f" elementwise.

julia> map((x) -> x * 2, [1, 2, 3])
3-element Array{Int64,1}:
2
4
6

julia> map(+, [1, 2, 3], [10, 20, 30])
3-element Array{Int64,1}:
11
22
33
In [142]:
lu(M2)
Out[142]:
(
100x100 Array{BigFloat,2}:
1e+00                                                                                 …  0e+00
1.838523316026038069423294824104585502107622762920382948468461321049017107287582e-01     0e+00
2.226652672640778428153578974542223057899272696712333628373710024305295178288089e-01     0e+00
8.800271373546953526824682655928247523713656239133309029931077636194740548240004e-01     0e+00
9.969711511643861050008990963774833601635694233722855443788571591033609498720476e-01     0e+00
8.979138730555442894457974595685742086077359813056139849438793026237408391767279e-01  …  0e+00
6.807588600102166833365124716208762430663471006732089118695641612340596096693672e-01     0e+00
2.455829705129951553304067727592573195211437393497849255275277512239547428646326e-01     0e+00
7.733884474442739333431587703752807939448278524827717594937087861232859025738641e-01     0e+00
6.369478597346308029171973720850203624729145464284684008762430813593740043172259e-01     0e+00
2.827829745702584616783062241061194670528384475163341040075224413128537647704774e-02  …  0e+00
6.8916911836155407719497226896317385232855759890401728866826844091322063237215e-02       0e+00
9.914847445449848621627477787268697917524046066273956341540063379389126840704449e-01     0e+00
⋮                                                                                     ⋱
5.917560795509736662262312292818648981295636865290101585946113867783160165303854e-01     0e+00
2.730343465699633349451394482373398269247854760360330831523299677598397581623309e-01     0e+00
9.039945797670267226102309624603231261343146722095674272138480977753157246173324e-01  …  0e+00
9.462017626922501325170145941450688057814574666304378855494052732254706598746676e-01     0e+00
7.922089379792833548468864680973299515938762464750620693002149619341420960071084e-01     0e+00
8.91312804253419886178066534318762609663940375583837777475766448455788684217046e-01      0e+00
2.792642705237073137742461219913623698283706391737491380867611670762794730952418e-05     0e+00
7.573650294135219377704801284400430317640792346704970287945898791774355476300744e-01  …  0e+00
4.340348981096360281012666931631459479552193700232478278658625428081540612698508e-01     0e+00
8.425932003310413942243844696491498599262760249784090231673454808048725464374257e-01     0e+00
9.574420989713105503105950862187424587135981099834153833267817582826759295038591e-01     0e+00
9.937128828323402433568347457741518522388968395097111212977949553756464564214356e-01     1e+00,

100x100 Array{BigFloat,2}:
9.994851507686071112601666754926554858684539794921875e-01  …   5.00779152589674314555168166407383978366851806640625e-02
0e+00                                                          5.957614820766671895663417312682526372739723842260882343182715867242722586063138e-01
0e+00                                                          4.469478342275836995433682687239315802830255734475567159873597681467845826134468e-01
0e+00                                                          1.594225612955579740227438886545435994704893579284040834538915577684978835022275e-01
0e+00                                                          3.51418102474667774543916344035850933636933755274101969213535739730486430086481e-01
0e+00                                                      …   3.853995991912961375481683208700128985955666356716160345936599990779107394373079e-01
0e+00                                                         -1.388928850457794922700969846901469586577590670720174042695508511464725140648597e-01
0e+00                                                         -1.963736253744149588710187786364289757425844420781275911168560619527357939135026e-01
0e+00                                                          3.178968402329430072315582253427327119635802392272492696354540300940096734477185e-01
0e+00                                                          2.176214988150537681248929128626937855865936488765697272966015579849791295536862e-01
0e+00                                                      …   5.119997557073573321754864325074991553008298165056610101108332375190160965916915e-01
0e+00                                                         -8.199420519824984131304944439473283677509641603892215267311260421539949314147843e-02
0e+00                                                         -1.100757926861985762949819942253220494124734219500225416371630465493477225668268e+00
⋮                                                          ⋱
0e+00                                                          3.119374583443352091823412074951378032702116370133229255033911885197311498449721e-01
0e+00                                                         -3.963473902224899280333406215892276084125185726495321949218801950999434601486697e-01
0e+00                                                      …   1.460022051064841242988743472828596892824793640959853693084700574539034836575632e+00
0e+00                                                          9.460562261692564028602511319807472584567333376119182633775041922970959644819617e-01
0e+00                                                          1.916100769270754349291340364996423402303124329778281922547647681942010509253265e+00
0e+00                                                          2.453777141130431399782611832707715011759461721288957685494547893432872020437905e+00
0e+00                                                          3.092241770787052786375610695487741355942081214332538992688100491955019679622272e-01
0e+00                                                      …  -7.690085518435571170493252506836527889003361725238261885731324700815770035046439e-01
0e+00                                                          7.057077397497909132059254043595925169793816268455912227899415362222622859241975e-01
0e+00                                                          2.11408055146617577269295905987997866661481241577267775907335892526116912022093e-01
0e+00                                                          4.72135706225358566629488251105981838448104777282753072505691201455089121862885e-01
0e+00                                                         -1.411605504403213576549449087599331540683320578509496538094218931547453928045177e+00,

[41,78,1,52,36,63,62,91,21,95  …  67,34,76,69,31,77,28,82,97,58])
In [143]:
lu(M)
Out[143]:
(
100x100 Array{Float64,2}:
1.0         0.0         0.0         0.0        …   0.0       0.0       0.0
0.183852    1.0         0.0         0.0            0.0       0.0       0.0
0.222665    0.333612    1.0         0.0            0.0       0.0       0.0
0.880027    0.908157    0.485849    1.0            0.0       0.0       0.0
0.996971    0.0121974   0.880719    0.482328       0.0       0.0       0.0
0.897914    0.688716   -0.173968   -0.535249   …   0.0       0.0       0.0
0.680759    0.256184    0.805169    0.0615926      0.0       0.0       0.0
0.245583    0.122523    0.852098    0.230656       0.0       0.0       0.0
0.773388    0.783385    0.294589    0.813216       0.0       0.0       0.0
0.636948    0.621283   -0.464593    0.0427875      0.0       0.0       0.0
0.0282783   0.895956   -0.28519    -0.723122   …   0.0       0.0       0.0
0.0689169   0.529095    0.623271    0.567317       0.0       0.0       0.0
0.991485    0.723194   -0.0744293   0.0926657      0.0       0.0       0.0
⋮                                              ⋱
0.591756    0.553754    0.585172    0.792699       0.0       0.0       0.0
0.273034    0.946579    0.312432    0.252827       0.0       0.0       0.0
0.903995    0.190557    0.113748   -0.253964   …   0.0       0.0       0.0
0.946202    0.644426    0.227142    0.306999       0.0       0.0       0.0
0.792209    0.616671    0.267607    0.875396       0.0       0.0       0.0
0.891313    0.895161    0.31949     0.0199976      0.0       0.0       0.0
2.79264e-5  0.187737    0.926009    0.460432       0.0       0.0       0.0
0.757365    0.538407    0.300644    0.140578   …   0.0       0.0       0.0
0.434035    0.293604    0.46559     0.452693       0.0       0.0       0.0
0.842593    0.325906    0.090703    0.620878       1.0       0.0       0.0
0.957442    0.628242    0.0057032  -0.198115      -0.163927  1.0       0.0
0.993713    0.465608   -0.431142   -0.125049      -0.617413  0.640633  1.0,

100x100 Array{Float64,2}:
0.999485  0.0617704  0.296413  …   0.536738     0.579976    0.0500779
0.0       0.964559   0.333478     -0.0592369   -0.0351511   0.595761
0.0       0.0        0.756701      0.751446     0.616732    0.446948
0.0       0.0        0.0          -0.572233    -0.401953    0.159423
0.0       0.0        0.0          -0.097152    -0.229394    0.351418
0.0       0.0        0.0       …   0.171652     0.103321    0.3854
0.0       0.0        0.0          -0.823456    -0.31806    -0.138893
0.0       0.0        0.0           0.346634     0.117767   -0.196374
0.0       0.0        0.0           0.114205     0.694677    0.317897
0.0       0.0        0.0           0.18621      0.454778    0.217621
0.0       0.0        0.0       …   0.875338     1.06396     0.512
0.0       0.0        0.0          -0.666298    -0.761167   -0.0819942
0.0       0.0        0.0          -0.463261    -0.622826   -1.10076
⋮                              ⋱
0.0       0.0        0.0           1.10803     -2.71365     0.311937
0.0       0.0        0.0           1.55302     -4.63257    -0.396347
0.0       0.0        0.0       …   0.00811131  -1.63731     1.46002
0.0       0.0        0.0           0.959823    -1.79163     0.946056
0.0       0.0        0.0           2.08064      0.0531956   1.9161
0.0       0.0        0.0           0.803113    -1.73451     2.45378
0.0       0.0        0.0          -1.11381      0.340109    0.309224
0.0       0.0        0.0       …   1.40266      0.277447   -0.769009
0.0       0.0        0.0          -1.36105      0.370272    0.705708
0.0       0.0        0.0          -1.64671      0.297193    0.211408
0.0       0.0        0.0           0.0          1.00489     0.472136
0.0       0.0        0.0           0.0          0.0        -1.41161  ,

[41,78,1,52,36,63,62,91,21,95  …  67,34,76,69,31,77,28,82,97,58])
In [144]:
?lu
Base.lu(A) -> L, U, p

Compute the LU factorization of "A", such that "A[p,:] = L*U".
In [145]:
methods(lu)
Out[145]:
2 methods for generic function lu:
In [146]:
# @edit lu(M)
In [147]:
?edit
Base.edit(file::String[, line])

Edit a file optionally providing a line number to edit at. Returns
to the julia prompt when you quit the editor.

Base.edit(function[, types])

Edit the definition of a function, optionally specifying a tuple of
types to indicate which method to edit.

## Interacting with the system¶

### Command-line arguments¶

A Julia script, similar to a Python script, is a sequence of Julia commands placed in a file, with the termination .jl.

From the command line, a script script.jl can be run as

julia script.jl arg1 arg2

where arg1 and arg2 are command-line arguments.

These command-line arguments to Julia scripts are placed in the variable ARGS as an array of strings.

### Files¶

Simple file input and output is easy:

In [148]:
outfile = open("test.txt", "w")
Out[148]:
IOStream(<file test.txt>)
In [149]:
for i in 1:10
println(outfile, "The value of i is $i") end close(outfile) In [150]: ;cat test.txt The value of i is 1 The value of i is 2 The value of i is 3 The value of i is 4 The value of i is 5 The value of i is 6 The value of i is 7 The value of i is 8 The value of i is 9 The value of i is 10 In [151]: infile = open("test.txt", "r") Out[151]: IOStream(<file test.txt>) In [152]: lines = readlines(infile) Out[152]: 10-element Array{Union(UTF8String,ASCIIString),1}: "The value of i is 1\n" "The value of i is 2\n" "The value of i is 3\n" "The value of i is 4\n" "The value of i is 5\n" "The value of i is 6\n" "The value of i is 7\n" "The value of i is 8\n" "The value of i is 9\n" "The value of i is 10\n" In [153]: map(split, lines) Out[153]: 10-element Array{Array{SubString{ASCIIString},1},1}: SubString{ASCIIString}["The","value","of","i","is","1"] SubString{ASCIIString}["The","value","of","i","is","2"] SubString{ASCIIString}["The","value","of","i","is","3"] SubString{ASCIIString}["The","value","of","i","is","4"] SubString{ASCIIString}["The","value","of","i","is","5"] SubString{ASCIIString}["The","value","of","i","is","6"] SubString{ASCIIString}["The","value","of","i","is","7"] SubString{ASCIIString}["The","value","of","i","is","8"] SubString{ASCIIString}["The","value","of","i","is","9"] SubString{ASCIIString}["The","value","of","i","is","10"] In [154]: [float(line[6]) for line in map(split, lines)] Out[154]: 10-element Array{Any,1}: 1.0 2.0 3.0 4.0 5.0 6.0 7.0 8.0 9.0 10.0 In [155]: x = rand(5,5) Out[155]: 5x5 Array{Float64,2}: 0.494374 0.548784 0.265116 0.14175 0.225112 0.610879 0.863226 0.90218 0.905241 0.0066329 0.352839 0.163394 0.514312 0.46365 0.750558 0.72891 0.539688 0.37556 0.0172088 0.0674222 0.069424 0.00747654 0.951597 0.984837 0.350702 In [156]: writedlm("random.txt", x) In [157]: ;cat random.txt .49437356567549995 .5487835147411582 .2651160558475778 .1417496876038895 .22511205695160474 .6108787919798604 .863226095224578 .9021796958394908 .9052411849611 .006632897854429665 .3528390706690816 .16339401668213926 .5143118341066373 .4636496583658838 .7505578231005474 .7289100209460473 .5396877767037889 .37555992341554356 .017208754116336156 .06742216219992869 .069423968121745 .007476542969398592 .9515968215572792 .9848367987309561 .35070195619498623 In [158]: y = readdlm("random.txt") # note that tab completion works for files Out[158]: 5x5 Array{Float64,2}: 0.494374 0.548784 0.265116 0.14175 0.225112 0.610879 0.863226 0.90218 0.905241 0.0066329 0.352839 0.163394 0.514312 0.46365 0.750558 0.72891 0.539688 0.37556 0.0172088 0.0674222 0.069424 0.00747654 0.951597 0.984837 0.350702 In [159]: ?writedlm Base.writedlm(f, A, delim='t') Write "A" (either an array type or an iterable collection of iterable rows) as text to "f" (either a filename string or an "IO" stream) using the given delimeter "delim" (which defaults to tab, but can be any printable Julia object, typically a "Char" or "String"). For example, two vectors "x" and "y" of the same length can be written as two columns of tab-delimited text to "f" by either "writedlm(f, [x y])" or by "writedlm(f, zip(x, y))". In [160]: ?readdlm Base.readdlm(source, delim::Char, T::Type, eol::Char; header=false, skipstart=0, skipblanks=true, use_mmap, ignore_invalid_chars=false, quotes=true, dims, comments=true, comment_char='#') Read a matrix from the source where each line (separated by "eol") gives one row, with elements separated by the given delimeter. The source can be a text file, stream or byte array. Memory mapped files can be used by passing the byte array representation of the mapped segment as source. If "T" is a numeric type, the result is an array of that type, with any non-numeric elements as "NaN" for floating-point types, or zero. Other useful values of "T" include "ASCIIString", "String", and "Any". If "header" is "true", the first row of data will be read as header and the tuple "(data_cells, header_cells)" is returned instead of only "data_cells". Specifying "skipstart" will ignore the corresponding number of initial lines from the input. If "skipblanks" is "true", blank lines in the input will be ignored. If "use_mmap" is "true", the file specified by "source" is memory mapped for potential speedups. Default is "true" except on Windows. On Windows, you may want to specify "true" if the file is large, and is only read once and not written to. If "ignore_invalid_chars" is "true", bytes in "source" with invalid character encoding will be ignored. Otherwise an error is thrown indicating the offending character position. If "quotes" is "true", column enclosed within double-quote () characters are allowed to contain new lines and column delimiters. Double-quote characters within a quoted field must be escaped with another double-quote. Specifying "dims" as a tuple of the expected rows and columns (including header, if any) may speed up reading of large files. If "comments" is "true", lines beginning with "comment_char" and text following "comment_char" in any line are ignored. Base.readdlm(source, delim::Char, eol::Char; options...) If all data is numeric, the result will be a numeric array. If some elements cannot be parsed as numbers, a cell array of numbers and strings is returned. Base.readdlm(source, delim::Char, T::Type; options...) The end of line delimiter is taken as "\n". Base.readdlm(source, delim::Char; options...) The end of line delimiter is taken as "\n". If all data is numeric, the result will be a numeric array. If some elements cannot be parsed as numbers, a cell array of numbers and strings is returned. Base.readdlm(source, T::Type; options...) The columns are assumed to be separated by one or more whitespaces. The end of line delimiter is taken as "\n". Base.readdlm(source; options...) The columns are assumed to be separated by one or more whitespaces. The end of line delimiter is taken as "\n". If all data is numeric, the result will be a numeric array. If some elements cannot be parsed as numbers, a cell array of numbers and strings is returned. #### Shelling out¶ In [161]: ;ls Example of defining a new type.ipynb fractal_example helllo_julia.png Index.ipynb Install IPython Notebook.md Install Julia.md Install Julia Packages.md Interoperability with Python and C.ipynb Introduction to Julia.ipynb intro.txt julia_tutorial.md Metaprogramming.ipynb Packages and graphics in Julia.ipynb Performance in Julia.ipynb Plan of the tutorial.ipynb random.txt README.md stuff.pdf test.txt Vector2d.jl What's next.ipynb Why Julia.ipynb Why Julia.slides.html In [162]: run(echo Hello) Hello ### Scientific computing¶ Linear algebra, FFT, random numbers, special functions. Packages for optimization, ODEs etc. ## Functions¶ Functions may be defined using the short syntax f(x) = 3x + 1 or using a longer form: In [163]: dup(x) = 2x Out[163]: dup (generic function with 1 method) In [164]: function duplicate(x) 2x # no explicit "return" needed end Out[164]: duplicate (generic function with 1 method) The last value computed in the function is automatically returned; no explicit return statement is required. In [165]: duplicate(x) = x^2 Out[165]: duplicate (generic function with 1 method) Every operator in Julia is a function. Functions are implemented by specifying their action on different types. Until now, we have written only functions that are generic, in the sense that they do not specify which type they accept, and as in Python they will work as long as the operations performed in them make sense for the input value: In [166]: duplicate(3), duplicate(3.5), duplicate(1+3im) Out[166]: (9,12.25,-8 + 6im) In [167]: duplicate("Hola") Out[167]: "HolaHola" In [168]: 2 * "Hola" * has no method matching *(::Int64, ::ASCIIString) while loading In[168], in expression starting on line 1 Note that string concatenation uses the * operator in Julia, instead of the + operator as in Python. Repeating a string is thus done by raising to an integer power: In [169]: "Hello"^2 Out[169]: "HelloHello" As a simple example, suppose that we wish to concatenate two strings. In Python we would write: In [170]: s1 = "Hello, " s2 = "David" s1 + s2 + has no method matching +(::ASCIIString, ::ASCIIString) while loading In[170], in expression starting on line 4 However, we see that in Julia, summation is not defined for strings. What is it defined for? In [171]: + Out[171]: + (generic function with 148 methods) We see that + is treated as a function, and that it has a multitude of methods, which, in Julia, are specialised versions of the function that act on different types: If we were unaware of the * operator for string concatenation, we could just define our own + for the concatenation of two strings: In [173]: +(s1::String, s2::String) = string(s1, s2) Out[173]: + (generic function with 149 methods) In [174]: "First" + " second" Out[174]: "First second" However, we cannot add a number to a string, since we have not (yet) defined it: In [175]: "The value of x is " + 3 + has no method matching +(::ASCIIString, ::Int64) while loading In[175], in expression starting on line 1 This we can also define, using the previous new definition: In [176]: +(s::String, x::Number) = s + "$(2x)"
Out[176]:
+ (generic function with 150 methods)
In [177]:
"The value of x is " + 3
Out[177]:
"The value of x is 6"
In [178]:
x = 3.5
"The value of x is " + x
Out[178]:
"The value of x is 7.0"
In [179]:
3 + "hello"
+ has no method matching +(::Int64, ::ASCIIString)

In fact, we can define the summation of a string with any other object:

In [180]:
+(s::String, x) = s + string(x)
Out[180]:
+ (generic function with 151 methods)
In [181]:
"Complex " + [3,4,5]
Out[181]:
"Complex [3,4,5]"
In [182]:
"a" + 3
Out[182]:
"a6"
In [183]:
Number
Out[183]:
Number
In [184]:
typeof(Number)
Out[184]:
DataType
In [185]:
super(Int64)
Out[185]:
Signed
In [186]:
super(Signed)
Out[186]:
Integer
In [187]:
super(Integer)
Out[187]:
Real
In [188]:
super(Real)
Out[188]:
Number
In [189]:
super(Number)
Out[189]:
Any

In this way, the concept of "function" is replaced by a "patchwork" of different definitions for objects of different types, easily modifiable by the user. This is also exactly the way to define "operator overloading" for user-defined types.

In the above, we also begin to see the power of multiple dispatch: we defined two methods (versions) of the function +, both with the same number but different types of arguments.

### User-defined types¶

A user-defined "composite type" is a collection of data. Unlike in Python, types do not "own" methods (functions internal to the type).

Rather, methods are defined separately, and are characterised by the types of all of their arguments; this is known as multiple dispatch. (Dispatch is the process of choosing which "version" of a given function to execute.)

A simple, but useful, example, is that of defining a 2D vector type. (See also the ImmutableArrays.jl package; fixed-size arrays will later be incorporated into base Julia.)

In [190]:
@which 3//4
Out[190]:
//(n::Integer,d::Integer) at rational.jl:15
In [191]:
Rational(3)
Out[191]:
3//1
In [192]:
6//4
Out[192]:
3//2
In [193]:
im*im
Out[193]:
-1 + 0im
In [194]:
immutable Vector2D   # type
x::Float64
y::Float64
end
In [195]:
v = Vector2D(3, 4)
w = Vector2D(5, 6)
Out[195]:
Vector2D(5.0,6.0)
In [196]:
v + w
+ has no method matching +(::Vector2D, ::Vector2D)
In [197]:
+(v::Vector2D, w::Vector2D) = Vector2D(v.x+w.x, v.y+w.y)
Out[197]:
+ (generic function with 152 methods)
In [198]:
v + w
Out[198]:
Vector2D(8.0,10.0)
In [199]:
*(v::Vector2D, α::Number) = Vector2D(v.x*α, v.y*α)
*(α::Number, v::Vector2D) = Vector2D(v.x*α, v.y*α)
Out[199]:
* (generic function with 126 methods)
In [200]:
v * 3.5
Out[200]:
Vector2D(10.5,14.0)
In [201]:
3.5 * v
Out[201]:
Vector2D(10.5,14.0)

Exercise: Define mathematical operations on Vector2D. Define a particle with position and velocity in 2D. Define function move that acts on a particle to move it over a time $\delta t$.

In [202]:
immutable VectorT   # type
p::Vector2D
v::Vector2D
end
In [203]:
a = VectorT(Vector2D(0,0),Vector2D(2,3))
Out[203]:
VectorT(Vector2D(0.0,0.0),Vector2D(2.0,3.0))
In [204]:
move(a::VectorT, t::Number) = a.p+a.v*t
Out[204]:
move (generic function with 1 method)
In [205]:
println("new position of particle $a with origin position$(a.p) and velocity $(a.v) at time 3 is$(move(a, 3))")
new position of particle VectorT(Vector2D(0.0,0.0),Vector2D(2.0,3.0)) with origin position Vector2D(0.0,0.0) and velocity Vector2D(2.0,3.0) at time 3 is Vector2D(6.0,9.0)

Here, we have used immutable instead of type for efficiency: the object is stored in an efficient packed form.

The equivalent of the Python __repr__ method for an object is to extend the show method:

In [206]:
?show
Base.show(x)

Write an informative text representation of a value to the current
output stream. New types should overload "show(io, x)" where the
first argument is a stream. The representation used by "show"
generally includes Julia-specific formatting and type information.
In [207]:
import Base.show

show(io::IO, v::Vector2D) = print(io, "[$(v.x),$(v.y)]")
Out[207]:
show (generic function with 94 methods)
In [208]:
v
Out[208]:
[3.0, 4.0]
In [209]:
+(v1::Vector2D, v2::Vector2D) = Vector2D(v1.x+v2.x, v1.y+v2.y)

*(v::Vector2D, lamb::Number)  = Vector2D(lamb*v.x, lamb*v.y)
Out[209]:
* (generic function with 126 methods)

We can confirm that the new method for the function + has indeed been defined:

In [211]:
x = Vector2D(3, 4)
y = Vector2D(5, 6)

x + y
Out[211]:
[8.0, 10.0]

### Parametrised types¶

Types may have a parameter, for example:

In [212]:
immutable Vector2D{T <: Real} #must same type T
x::T
y::T
end
invalid redefinition of constant Vector2D

T is a type parameter. The expression T <: Real means that T must be a subtype of the abstract type Real. We can investigate the type hierarchy with the super function:

In [213]:
Integer
Out[213]:
Integer
In [214]:
super(Integer)
Out[214]:
Real
In [215]:
super(Real)
Out[215]:
Number
In [216]:
super(Number)
Out[216]:
Any
In [217]:
v = Vector2D(3., 4.)
Out[217]:
[3.0, 4.0]
In [218]:
v = Vector2D(3, 4.)
Out[218]:
[3.0, 4.0]

Here, the types of the two arguments were different, so there is no match for the type signature.

In [219]:
v = Vector2D(3//4, 5//6)
Out[219]:
[0.75, 0.8333333333333334]
In [220]:
import Base.show
show{T}(io::IO, v::Vector2D{T}) = print(io, "[$(v.x),$(v.y)]")
too many parameters for type Vector2D
In [221]:
v
Out[221]:
[0.75, 0.8333333333333334]

We can define outer constructors, defined outside the type definition itself, which allow other ways of constructing the object:

In [222]:
Vector2D{T}(x::T) = Vector2D(x, x)
Out[222]:
Vector2D (constructor with 3 methods)
In [223]:
Vector2D(3)
Out[223]:
[3.0, 3.0]

### Example: a simple type for a collection of particles¶

Let's define a particle:

In [224]:
type Particle
position::Vector2D{Float64}
velocity::Vector2D{Float64}
end
move(p::Particle, dt::Real) = p.position += p.velocity * dt
too many parameters for type Vector2D
In [225]:
show(io::IO, p::Particle) = print(io, "pos: $(p.position); vel:$(p.velocity)")
Particle not defined
In [226]:
p = Particle(Vector2D(0.,0.), Vector2D(1.,1.))
+{T}(v1::Vector2D{T}, v2::Vector2D{T}) = Vector2D{T}(v1.x+v2.x, v1.y+v2.y)

*{T}(v::Vector2D{T}, lamb::Number)  = Vector2D{T}(lamb*v.x, lamb*v.y)
Particle not defined
In [227]:
move(p, 0.1)
p not defined
In [228]:
p
p not defined

Now we can define a gas as a collection of particle:

In [229]:
Int
Out[229]:
Int64
In [230]:
type Gas
particles::Vector{Particle}  # Array{Particle, 1}

function Gas(N::Integer)
parts = [Particle(Vector2D(rand(2)...), Vector2D(rand(2)...)) for i in 1:N]
new(parts)
end
end
Particle not defined
In [231]:
show(io::IO, g::Gas) = for i in 1:length(g.particles); \
println(io, "Particle $i:$(g.particles[i])"); end
Gas not defined
In [232]:
g = Gas(10)
Gas not defined
In [233]:
move(g, 1)
g
move has no method matching move(::Function, ::Int64)
In [234]:
function move(g::Gas, dt::Number)
for particle in g.particles
move(particle, dt)
end
end
Gas not defined
In [235]:
move(g, 1)
move has no method matching move(::Function, ::Int64)
In [236]:
g
Out[236]:
g (generic function with 1 method)

### Matrix multplication¶

In [237]:
M = [1 1; 0 1]
Out[237]:
2x2 Array{Int64,2}:
1  1
0  1
In [238]:
M * M
Out[238]:
2x2 Array{Int64,2}:
1  2
0  1
In [239]:
M^10
Out[239]:
2x2 Array{Int64,2}:
1  10
0   1
In [240]:
M^(-1)
Out[240]:
2x2 Array{Float64,2}:
1.0  -1.0
0.0   1.0
In [241]:
inv(M)
Out[241]:
2x2 Array{Float64,2}:
1.0  -1.0
0.0   1.0
In [242]:
?inv
Base.inv(M)

Matrix inverse
In [243]:
det(M)
Out[243]:
1
In [244]:
norm(M)
Out[244]:
1.618033988749895
In [245]:
?norm
Base.norm(A[, p])

Compute the "p"-norm of a vector or the operator norm of a matrix
"A", defaulting to the "p=2"-norm.

For vectors, "p" can assume any numeric value (even though not
all values produce a mathematically valid vector norm). In
particular, "norm(A, Inf)" returns the largest value in
"abs(A)", whereas "norm(A, -Inf)" returns the smallest.

For matrices, valid values of "p" are "1", "2", or "Inf".
(Note that for sparse matrices, "p=2" is currently not
implemented.) Use "vecnorm()" to compute the Frobenius norm.
In [246]:
lamb, vv = eig(M)
Out[246]:
([1.0,1.0],
2x2 Array{Float64,2}:
1.0  -1.0
0.0   2.22045e-16)
In [247]:
Symmetric(M)
Out[247]:
2x2 Symmetric{Int64,Array{Int64,2}}:
1  1
1  1
In [248]:
?Symmetric
DataType   : Symmetric{T,S<:AbstractArray{T,2}} (constructor with 2 methods)
supertype: AbstractArray{T,2}
fields   : (:data,:uplo)
In [249]:
methods(Symmetric)
Out[249]:
2 methods for generic function Symmetric:

### Eigenvalues of a random matrix¶

Exercise: Generate a random matrix of $1000 \times 1000$ gaussian random variates. Calculate its eigenvalues and plot them. Calculate the differences between consecutive eigenvalues and plot a histogram of them.

In [250]:
N = 1000
M = randn(N, N)
M = Symmetric(M);
In [251]:
typeof(M)
Out[251]:
Symmetric{Float64,Array{Float64,2}} (constructor with 1 method)
In [252]:
@which dot([3], [4])
Out[252]:
dot(x::AbstractArray{T,1},y::AbstractArray{T,1}) at linalg/matmul.jl:49
In [253]:
help("dot")
Base.LinAlg.BLAS.dot(n, X, incx, Y, incy)

Dot product of two vectors consisting of "n" elements of array
"X" with stride "incx" and "n" elements of array "Y" with
stride "incy".
In [254]:
@time lamb, vv = eig(M);
elapsed time: 0.730055745 seconds (24964608 bytes allocated, 6.74% gc time)
In [255]:
using PyPlot
In [256]:
plot(lamb, "o-")
Out[256]:
1-element Array{Any,1}:
PyObject <matplotlib.lines.Line2D object at 0x7ff8f8bd3390>
In [257]:
PyPlot.figure(figsize=(8,4))
plot(lamb, "o-")
Out[257]:
1-element Array{Any,1}:
PyObject <matplotlib.lines.Line2D object at 0x7ff8f8a19310>
In [258]:
differences = diff(lamb);
In [259]:
plot(differences, "o-")
Out[259]:
1-element Array{Any,1}:
PyObject <matplotlib.lines.Line2D object at 0x7ff8f89d3650>
In [260]:
h = hist(differences, 100)
Out[260]:
(0.0:0.01:0.93,[8,16,23,48,47,54,59,56,54,66  …  0,0,0,0,0,0,0,0,0,1])
In [261]:
plot(collect(h[1][1:end-1]), h[2])
Out[261]:
1-element Array{Any,1}:
PyObject <matplotlib.lines.Line2D object at 0x7ff8f88ed850>

### requires vs include vs using vs import¶

In [262]:
g.particles
type Function has no field particles