User Tools

Site Tools


arraysmatrix

Differences

This shows you the differences between two versions of the page.

Link to this comparison view

Both sides previous revision Previous revision
Next revision
Previous revision
arraysmatrix [2018/12/11 02:39]
julia [For Loops]
arraysmatrix [2018/12/28 12:00] (current)
Line 11: Line 11:
  
 ```juliarepl ```juliarepl
-julia> pkgchk( [ "​julia"​ => v"1.0.2" ] )+julia> pkgchk.( [ "​julia"​ => v"1.0.3" ] );
 ``` ```
  
Line 20: Line 20:
  
 A 2-D array can be created by specifying all elements in blank separated columns and semicolon (not comma!) separated rows: A 2-D array can be created by specifying all elements in blank separated columns and semicolon (not comma!) separated rows:
- 
-FIXME (Andreas:) Also mention that you can use newline instead of semicolon. It can be very convenient when copy-pasting. 
  
 ```juliarepl ```juliarepl
Line 29: Line 27:
  ​5 ​ 6  7  8  ​5 ​ 6  7  8
 ``` ```
 +
 +- A newline can be used instead of the semi-colon
  
 WARNING Watch the spaces: ​ `x= [ 1 -2 ]` is a 1x2 array (matrix). ​ In contrast, `x= [ 1-2 ]` is a 1-element vector. WARNING Watch the spaces: ​ `x= [ 1 -2 ]` is a 1x2 array (matrix). ​ In contrast, `x= [ 1-2 ]` is a 1-element vector.
Line 35: Line 35:
  
 ```juliarepl ```juliarepl
-julia> x= [ 1, 2, 3 ] ; y= [ 10, 11, 12 ]; z= [ x, y ]+julia> x= [ 1, 2, 3 ] ; y= [ 10, 11, 12 ]; 
 + 
 +julia> ​[ x, y ]
 2-element Array{Array{Int64,​1},​1}:​ 2-element Array{Array{Int64,​1},​1}:​
  [1, 2, 3]  [1, 2, 3]
Line 41: Line 43:
 ``` ```
  
-### hcat/​vcat/​hvcat+### Cat (hcat/​vcat/​hvcat) and Bracket Alternatives
  
 Vectors are created with comma-separated lists. ​ Matrices are created with space and semicolon separated lists. Vectors are created with comma-separated lists. ​ Matrices are created with space and semicolon separated lists.
  
-The parenthesis notation above is just shortcuts.+The parenthesis notation above are "just" ​shortcuts.
  
 ```juliarepl ```juliarepl
Line 54: Line 56:
 true true
  
-julia> hvcat((3,​3,​3),​ 1,​2,​3,​4,​50,​6,​7,​8,​9)+julia> hvcat( (3,3,3), 1,​2,​3,​4,​50,​6,​7,​8,​9 )
 3×3 Array{Int64,​2}:​ 3×3 Array{Int64,​2}:​
  ​1 ​ 2  3  ​1 ​ 2  3
Line 70: Line 72:
 ``` ```
  
-FIXME (Andreas:) I think you don't want to explicitly call `xcat` in most cases. Instead, you'd usually prefer the square brackets syntax for concatenation inspired by Matlab. 
  
  
Line 90: Line 91:
 ## Initializing a Matrix with Same Values ## Initializing a Matrix with Same Values
  
-`fill` can be used to initialize arbitrarily-dimensional arrays:+`fill()` can initialize arbitrarily-dimensional arrays:
  
 ```juliarepl ```juliarepl
-julia> fill( NaN , (3,2) )+julia> fill( NaN , (3,2) ) # or use fill( NaN, 3, 2 )
 3×2 Array{Float64,​2}:​ 3×2 Array{Float64,​2}:​
  ​NaN ​ NaN  ​NaN ​ NaN
Line 100: Line 101:
 ``` ```
  
-There are type-specialized versions for `zeros(Float64,​ 3, 3)`  and `ones(Bool, 3, 3)`, but not `fill()`.+There are type-specialized versions for `zeros(Float64,​ 3, 3)`  and `ones(Bool, 3, 3)`.
  
 You could also create *uninitialized* matrices. ​ Even if you immediately assign to your matrix afterwards, such uninitialized variables are rarely a good idea, and the cost of initialization with nada is trivial. ​ Chances are that sooner or later, you may forget. ​ Please always initialize your variables! You could also create *uninitialized* matrices. ​ Even if you immediately assign to your matrix afterwards, such uninitialized variables are rarely a good idea, and the cost of initialization with nada is trivial. ​ Chances are that sooner or later, you may forget. ​ Please always initialize your variables!
Line 107: Line 108:
 ### Diagonal Matrices ### Diagonal Matrices
  
-#### Arbitrary Diagonal+The identity matrix has special support: 
 + 
 +```juliarepl 
 +julia> using LinearAlgebra 
 + 
 +julia> Matrix( 2.0I, 3, 3 ) ## An '​I'​ is following '​2.0'​ creates a particular diagonal matrix 
 +3×3 Array{Float64,​2}:​ 
 + ​2.0 ​ 0.0  0.0 
 + ​0.0 ​ 2.0  0.0 
 + ​0.0 ​ 0.0  2.0 
 + 
 +julia> Matrix( I, 3, 3 ) ## I is shortcut for 1.0I 
 +3×3 Array{Float64,​2}:​ 
 + ​1.0 ​ 0.0  0.0 
 + ​0.0 ​ 1.0  0.0 
 + ​0.0 ​ 0.0  1.0 
 + 
 +julia> Matrix( 10.0I, 3, 3 ) + I ## I is an identity matrix that reshapes itself to whatever is needed 
 +3×3 Array{Float64,​2}:​ 
 + ​11.0 ​  ​0.0 ​  0.0 
 +  0.0  11.0   0.0 
 +  0.0   ​0.0 ​ 11.0 
 + 
 +``` 
 + 
 +For arbitrary diagonal matrices, you can use
  
 ```juliarepl ```juliarepl
Line 120: Line 146:
 ``` ```
  
-#### The Identity Matrix 
  
-```juliarepl 
-julia> using LinearAlgebra 
- 
-julia> Matrix( 1.0I, 3, 3 ) ## look closely. ​ An '​I'​ is following '​1.0'​ 
-3×3 Array{Float64,​2}:​ 
- ​1.0 ​ 0.0  0.0 
- ​0.0 ​ 1.0  0.0 
- ​0.0 ​ 0.0  1.0 
-``` 
- 
-FIXME! (Andreas) You might want to explain `I` since in most cases it is preferred over `Matrix( 1.0I, 3, 3 )`. 
  
 ## Pushing a Vector onto an Array ## Pushing a Vector onto an Array
 +
 +Remember: A matrix is not like a vector of vectors.
  
 ### Pushing Vectors onto Vectors ### Pushing Vectors onto Vectors
Line 158: Line 174:
 ``` ```
  
-### Appending Vectors onto Matrices +### Side-Append Column Vectors ​to Matrix
- +
-#### Side-Append Column Vectors+
  
 ```juliarepl ```juliarepl
Line 185: Line 199:
 ``` ```
  
-#### Bottom-Append Row Vectors+### Bottom-Append Row Vectors ​to Matrix
  
 ```juliarepl ```juliarepl
Line 227: Line 241:
  
  
-## Matrix Function From Two Vectors (for 3D or Contourplots)+## Matrix Function ​Over Grid From Two Vectors (for 3D or Contourplots)
  
-FIXME (Andreas:) You'd have to explain what is going on here.+Julia makes it easy to create a rectangular grid with values Here is how the grid will look like:
  
 ```juliarepl ```juliarepl
-julia> x= [1:5;]; y= [3:9;]; f(x,y)= x^2*sqrt(y);+julia> x= 1:5; y= 3:9; g(x,​y)= ​(x,y); ## show off the grid that we will generate
  
-julia> ​f.(x,​y'​) ​       ## note element-wise use with transpose+julia> ​g.(x,​y'​)         ## note element-wise use with transpose 
 +5×7 Array{Tuple{Int64,​Int64},​2}:​ 
 + (1, 3)  (1, 4)  (1, 5)  (1, 6)  (1, 7)  (1, 8)  (1, 9) 
 + (2, 3)  (2, 4)  (2, 5)  (2, 6)  (2, 7)  (2, 8)  (2, 9) 
 + (3, 3)  (3, 4)  (3, 5)  (3, 6)  (3, 7)  (3, 8)  (3, 9) 
 + (4, 3)  (4, 4)  (4, 5)  (4, 6)  (4, 7)  (4, 8)  (4, 9) 
 + (5, 3)  (5, 4)  (5, 5)  (5, 6)  (5, 7)  (5, 8)  (5, 9) 
 + 
 +``` 
 + 
 +The prime on the y defines the transpose (or adjoint). ​ Instead of a `g()` function that returns tuples to illustrate what will be used in each array cell, we can use a simple function: 
 + 
 +```juliarepl 
 +julia> x= 1:5; y= 3:9; f(x,y)= x^2*sqrt(y);​ ##​ apply a useful function 
 + 
 +julia> f.(x,​y'​)
 5×7 Array{Float64,​2}:​ 5×7 Array{Float64,​2}:​
   1.73205 ​  ​2.0 ​  ​2.23607 ​  ​2.44949 ​  ​2.64575 ​  ​2.82843 ​  3.0   1.73205 ​  ​2.0 ​  ​2.23607 ​  ​2.44949 ​  ​2.64575 ​  ​2.82843 ​  3.0
Line 250: Line 279:
  
 ```juliarepl ```juliarepl
-julia> M= [r * c for r in 1:3, c in 1:6]+julia> M= [r * c for r in 1:3, c in 1:6] ## a sample two-dimensional comprehension
 3×6 Array{Int64,​2}:​ 3×6 Array{Int64,​2}:​
  ​1 ​ 2  3   ​4 ​  ​5 ​  6  ​1 ​ 2  3   ​4 ​  ​5 ​  6
Line 256: Line 285:
  ​3 ​ 6  9  12  15  18  ​3 ​ 6  9  12  15  18
  
-julia> length(M) +julia> ndims(M) ## array as two dimensions
-18 +
- +
-julia> ndims(M)+
 2 2
  
-julia> size(M) ​   ## the number of rows (aka nrows, 3) and the number of columns (aka ncols, 6)+julia> size(M) ## number of rows (aka nrows, 3) and number of columns (aka ncols, 6)
 (3, 6) (3, 6)
  
-julia> size(M,2)+julia> size(M,2) ## the second dimension
 6 6
  
-julia> eachindex(M)+``` 
 + 
 + 
 +## Linear Indexing (Folded Matrix as Vector) 
 + 
 +Julia can interpret matrices as long vectors, placed in column-row order. ​ This is called *linear indexing*:​ 
 + 
 + 
 +```juliarepl 
 +julia> M= [r * c for r in 1:3, c in 1:6]; ## the same two-dimensional matrix 
 + 
 +julia>  
 + 
 + 
 + 
 +julia> length(M) ##​ the total number of elements 
 +18 
 + 
 +julia> ( M[5], M[8], M[11] ) ## 5=(2,2); 8=(2,3); 11=(2,4) 
 +( 4, 6, 8 ) 
 + 
 +julia> eachindex(M) ## An iterator (OneTo). (could be collect()ed)
 Base.OneTo(18) Base.OneTo(18)
  
-julia> ​collect(eachindex(M)) == collect([1:length(M);])+julia> eachindex(M) == 1:length(M)
 true true
  
-julia> ( stride(M,1) , stride(M,2) )+julia> ( stride(M,1) , stride(M,2) ) ## num indexes to advance in vector M to advance in matrix M direction 1 or 2
 (1, 3) (1, 3)
  
 ``` ```
- 
-FIXME (Andreas:) I think more explanation is needed here. What is `OneTo` and what is `stride`? 
  
  
-## Accessing Rows, Columns, or Cells of a Matrix+### Mixing Linear and Non-Linear Access
  
 ```juliarepl ```juliarepl
Line 290: Line 335:
  ​5 ​ 6  7  8  ​5 ​ 6  7  8
  
-julia> M[:]                ## matrix ​accessed ​as vector+julia> M[:]                ## matrix ​linear access ​as vector
 8-element Array{Int64,​1}:​ 8-element Array{Int64,​1}:​
  1  1
Line 301: Line 346:
  8  8
  
-julia> M[:,​1] ​             ## column.    ​(same result ​as x[ [1,2] ] in this case)+julia> M[:,​1] ​              ## column (same as x[ [1,2] ])
 2-element Array{Int64,​1}:​ 2-element Array{Int64,​1}:​
  1  1
  5  5
  
-julia> M[1,:​] ​             ## row+julia> M[1,:​] ​              ## row access
 4-element Array{Int64,​1}:​ 4-element Array{Int64,​1}:​
  1  1
Line 313: Line 358:
  4  4
  
-julia> M[2,​3] ​             ## cell+julia> M[2,​3] ​              ## cell access
 7 7
  
-julia> M[ 2, [3,4] ]       ​## arbitrary ​columns ​and rows+julia> M[ 2, [3,4] ]       ## arbitrary ​column ​and row access
 2-element Array{Int64,​1}:​ 2-element Array{Int64,​1}:​
  7  7
Line 322: Line 367:
  
  
-julia> M[ :, end-1] ​       ## `end` is useful+julia> M[ :, end-1] ​        ## `end` is often easier
 2-element Array{Int64,​1}:​ 2-element Array{Int64,​1}:​
  3  3
Line 333: Line 378:
 ``` ```
  
-To obtain an expression that you can paste into the REPL, use `@show M`. 
  
  
-## Direct Access to Matrix as Folded Vector+## REPL-Pastable Expression
  
 ```juliarepl ```juliarepl
-julia> M= [ 1 4 ; 5 7 8 ]; M[2]       ## internally stored in column-row order +julia> M= [r * c for r in 1:3, c in 1:6];        ## the two-dimensional matrix 
-5+ 
 +julia> @show(M); 
 +M = [1 2 3 4 6; 2 4 6 8 10 12; 3 6 9 12 15 18] 
 ``` ```
  
-FIXME (Andreas:) Mention that this is usually called *linear indexing* in Julia.+ 
  
 ## Converting an X-Dimensional Array into a Different (e.g., Y-Dimensional) Array ## Converting an X-Dimensional Array into a Different (e.g., Y-Dimensional) Array
Line 361: Line 409:
  
  
-#### Warning: Reshape is an Alias 
  
-If you change a value in the reshaped array, the original array will still be affected, too.  If this is not desired, work with a `copy(v)`.+WARNING **Reshape is an Alias**: ​If you change a value in the reshaped array, the original array will still be affected, too.  If this is not desired, work with a `copy(v)`.
  
 ```juliarepl ```juliarepl
-julia> M= [1,​2,​3,​4,​5,​6]; ​ b= reshape(M, 2, 3)+julia> M= [1,​2,​3,​4,​5,​6];​ ## the original matrix M 
 + 
 +julia> ​b= reshape(M, 2, 3) ## b is only an alias to M
 2×3 Array{Int64,​2}:​ 2×3 Array{Int64,​2}:​
  ​1 ​ 3  5  ​1 ​ 3  5
  ​2 ​ 4  6  ​2 ​ 4  6
  
-julia> b[1]= 99; M+julia> b[1]= 99; ## assigning to b changes ​M, too! 
 + 
 +julia> M ## as you can see here
 6-element Array{Int64,​1}:​ 6-element Array{Int64,​1}:​
  99  99
Line 408: Line 459:
  
 ``` ```
- 
- 
  
  
Line 435: Line 484:
 ``` ```
  
-FIXME (Andreas:) Mentioned that often it is not necessary to materialize the grid with `collect` ​and that it can save a lot of memory.+- It may not be necessary to materialize the grid with `collect`, which can save memory. 
  
 ## Learning the Dimension of Higher-Dimensional Arrays ## Learning the Dimension of Higher-Dimensional Arrays
Line 482: Line 532:
 ### For Loops ### For Loops
  
-The "​naive"​ and good-enough ​code is:+The standard ​code is:
  
 ```juliarepl ```juliarepl
Line 489: Line 539:
 julia> xset= 10:11; yset= 12:13; i= 1; julia> xset= 10:11; yset= 12:13; i= 1;
  
-julia> m= Matrix{Float64}(undef, length(xset)*length(yset),​ 3);+julia> m= fillNaN, length(xset)*length(yset),​ 3 )
  
 julia> i julia> i
Line 499: Line 549:
  m[i,:]= [ x, y, f(x,y) ]  m[i,:]= [ x, y, f(x,y) ]
  i+=1  i+=1
-     end#​for ​x +     end#​for ​y# 
- end#for##+ end#​for ​x##
  
 julia> m julia> m
Line 510: Line 560:
 ``` ```
  
-For speed, it can be improved to+Possible ​speed enhancements:​
  
-```julia +- use `Matrix{Float64}(undef, length(xset)*length(yset), ​3);to work with unitialized matrix
-    @inbounds for x in Float64.(xset) ​ +
-        @simd for y in Float64.(yset) +
-            m[:,i] .= ( x, y, f(x,y+
-            i += 1 +
-        end +
-    end +
-```+
  
-FIXME (Andreas:) I doubt that `@simdhelps here so maybe better to leave it out. Also, explain ​`@inbounds`. In particular what you risk to crash Julia when using it. +- place `@inboundsbefore the `forstatement tells Julia not to check the bounds.  ​Howeverthis risks crashing Julia and probably works only inside functions.
- +
-Note that `Float64.` ​to help SIMD.  ​The following is also surprisingly fast: +
- +
-```julia +
- m = Matrix{Float64}(undef3, length(xset)*length(yset)) +
- i = 1 +
- @inbounds for x in xset, y in yset  +
- m[1,i] = x +
- m[2,i] = y +
-        m[3,i] = f(x,y) +
- i += 1 +
- end +
-```+
  
 +- use `@simd @inbounds` before the `for` statement to ask for modestly advanced Intel CPU vector features.
  
  
Line 570: Line 601:
 using Base.Iterators using Base.Iterators
 reduce(hcat,​ [[x, y, f(x, y)] for (y, x) in Iterators.product(yset,​ xset)]) ## 8.531 μs reduce(hcat,​ [[x, y, f(x, y)] for (y, x) in Iterators.product(yset,​ xset)]) ## 8.531 μs
-``````+```
  
  
Line 578: Line 609:
 ## Broadcasting ## Broadcasting
  
-[[functions|dot_postfix_functions]] operate on scalars inside arrays. ​ To operate on arrays inside arrays (e.g., vectors inside matrices),+[[funother|dot_postfix_functions]] operate on scalars inside arrays. ​ To operate on arrays inside arrays (e.g., vectors inside matrices),
  
 ```juliarepl ```juliarepl
-julia> broadcast(/,​ [1 2 3 4 ; 5 6 7 8], [10,20])+julia> broadcast( / , [1 2 3 4 ; 5 6 7 8], [10,20] ) ## shorthand for [1 2 3 4 ; 5 6 7 8 ] ./ [10,20]
 2×4 Array{Float64,​2}:​ 2×4 Array{Float64,​2}:​
  ​0.1 ​  ​0.2 ​ 0.3   0.4  ​0.1 ​  ​0.2 ​ 0.3   0.4
  ​0.25 ​ 0.3  0.35  0.4  ​0.25 ​ 0.3  0.35  0.4
  
-julia> broadcast(/,​ [10,20], [1 2 3 4 ; 5 6 7 8])+julia> broadcast( / , [1 2 3 4 ; 5 6 7 8], [10,20] ) ==  [1 2 3 4 ; 5 6 7 8 ] ./ [10,20] 
 +true 
 + 
 +julia> broadcast( / , [10,20], [1 2 3 4 ; 5 6 7 8] ) ## shorthand for [10,20] ./ [1 2 3 4 ; 5 6 7 8 ]
 2×4 Array{Float64,​2}:​ 2×4 Array{Float64,​2}:​
  ​10.0 ​ 5.0      3.33333 ​ 2.5  ​10.0 ​ 5.0      3.33333 ​ 2.5
   4.0  3.33333 ​ 2.85714 ​ 2.5   4.0  3.33333 ​ 2.85714 ​ 2.5
-``` 
  
 +julia> broadcast( / , [10,20], [1 2 3 4 ; 5 6 7 8] ) == [10,20] ./ [1 2 3 4 ; 5 6 7 8 ]
 +true
 +
 +```
  
 +The dot function versions are faster than the broadcasts, because the dot functions know how to "​fuse"​ operations. ​ That is, they do not create temporary intermediate matrices---they fuse all the necessary operations into each destination cell.
  
  
  
-## Removing Sets of Rows Based on Condition (e.g., ​with NaN)+## Removing Sets of Rows Based on Condition (e.g., ​containing at least one NaN)
  
 ```juliarepl ```juliarepl
arraysmatrix.txt · Last modified: 2018/12/28 12:00 (external edit)