1 of 57

Programming Basics for Logistics Algorithms: Lecture 5

Tatiana Polishchuk

Associate Professor, 

Linköping University, KTS

www.itn.liu.se/~tatpo46/

1

2 of 57

Course Overview

  • Introduction to programming
  • Variables and data types
  • Algorithms and scripts
  • Input/output statements
  • Conditional statements
  • Loops, nested loops
  • Vectors
  • Matrices and graphs
  • Plotting: graphical functions
  • Good programming practices
  • Data management

2

3 of 57

Course Overview

  • Introduction to programming
  • Variables and data types
  • Algorithms and scripts
  • Input/output statements
  • Conditional statements
  • Loops, nested loops
  • Vectors
  • Matrices and graphs
  • Plotting: graphical functions
  • Good programming practices
  • Data management

3

4 of 57

Flashback: Lecture 4

4

Questions?

5 of 57

Plotting

5

x = 11;

y = 18;

plot(x,y,'r*')

xlabel('Time')

ylabel('Temperature')

title('Time and Temp')

axis([x-2 x+2 y-5 y+5])

6 of 57

Plot Options

6

  • Can change the line color, marker style, and line style by adding a string argument

  • Can plot without connecting the dots by omitting line style argument

  • Look at help plot for a full list of colors, markers, and line styles

7 of 57

Saving Figures

7

  • Figures can be saved in many formats. The common ones are:

.fig preserves all information

.bmp uncompressed image

.eps high-quality scaleable format

.pdf compressed image

8 of 57

3D Line Plots

8

  • We can plot in 3 dimensions just as easily as in 2D

» t=0:0.001:4*pi;

» x = sin(t);

» y = cos(t);

» plot3(x,y,t,'k','LineWidth',2);

  • Use tools on figure to rotate it

  • Can set limits on all 3 axes

» xlim, ylim, zlim

9 of 57

Specialized Plotting Functions

9

  • MATLAB has a lot of specialized plotting functions
  • polar-to make polar plots

» polar (0:0.01:2*pi, cos ((0:0.01:2*pi)*2));

  • bar-to make bar graphs

» bar (1:10, rand(1,10));

  • quiver-to add velocity vectors to a plot

» [X,Y] = meshgrid (1:10,1:10);

» quiver(X,Y, rand(10),rand(10));

  • stairs-plot piecewise constant functions

» stairs(1:10, rand(1,10));

  • fill-draws and fills a polygon with specified vertices

» fill ([0 1 0.5], [0 0 1], ‘r’);

  • see help on these functions for syntax
  • doc specgraph – for a complete list

10 of 57

Visualizing Graphs

digraph objects represent directed graphs, which have directional edges connecting the nodes

clear ;

A = [

0 1 0 1;

1 0 1 0;

0 1 0 1;

1 0 1 0

];

figure;

G = digraph (A);

plot (G);

10

11 of 57

Visualizing Graphs

Defining edge list - another way of creating graphs

E_list = [

1, 2;

1, 3;

2, 3

]; % edge list

E_weights = [1 1 2];

G = graph(E_list(:,1), E_list(:,2), E_weights);

% to add labels to the nodes:

G.Nodes.Name = { 'Olle ' 'Sara ' 'Karl ' }';

% to mark weights on the graph:

plot (G , 'EdgeLabel', G.Edges.Weight);

11

12 of 57

Next: Data Management

12

Questions?

13 of 57

13

File Operations

  • Read from the file: extract data from files

  • Write to the file: write data to the file starting from the very beginning

  • Append to the file: write data to the file, adding information to what is already there

  • Display content of the file:show file content

14 of 57

14

Read from the file

  • Read from the file: extract data from the file

load filename

  • assigns file content to the array filename
  • use files with extensions .dat, .txt, .mat, etc.
  • only with the data of the same format!
  • only data in matrix format (same number of elements in each row)

15 of 57

15

File Operations

  • Write to the file: write data to the file starting from the very beginning

save filename matrixname -ascii

  • creates file filename
  • uses file extensions .txt, .dat, ect.
  • .mat used to store variables from your workspace
  • writes the data from array matrixname
  • if file not empty -- overwrites the data!!
  • ‘-ascii’ provides correct data format

16 of 57

16

File Operations

  • Append to the file: write data to the file, adding information to what is already there

save filename matrixname -ascii -append

  • adds the data from array matrixname to the end of file filename

17 of 57

17

File Operations

  • Display content of the file: show file content

type filename

Example:

» type testfile.dat

2 4 6 1 0 0 4

2 3 4 1 2 3 4

1 1 2 3 4 5 1

1 2 5 6 2 3 5

18 of 57

18

Advanced File Operations

  • Read table of mixed data types

T = readtable(filename)

  • creates a table T by reading column oriented data from a file
  • determines the file format from the file extension:

.txt, .dat, or .csv for delimited text files

.xls, .xlsb, .xlsm, .xlsx, .xltm, .xltx, or .ods for spreadsheet files

  • creates one variable in T for each column in the file and reads variable names from the first row of the file
  • by default, creates variables that have data types that are appropriate for the data values detected in each column of the input file

19 of 57

19

Advanced File Operations

  • Read table of mixed data types

T = readtable(filename)

  • To access variables representing columns of the table:

T.columnname

  • Write of different data types to a file

writetable(T, filename)

20 of 57

20

Cell Arrays

  • Elements in cell arrays are containers that can store different types of values

  • stores pointers do different data
  • unique data structure for Matlab
  • use {} to create such arrays

Example:

» cellrowvec = {23, 'a', 1:2:9, 'hello'}

cellrowvec =

[23] 'a' [1x5 double] 'hello'

21 of 57

21

Cell Arrays: Indexing

  • Indexing is different from the data arrays

Example:

» cellrowvec = {23, 'a', 1:2:9, 'hello'}

cellrowvec =

[23] 'a' [1x5 double] 'hello'

» cellrowvec(3)

ans =

1×1 cell array

{1×5 double}

» cellrowvec{3}

ans =

1 3 5 7 9

Use curly braces {}

to access values stored in cell array!

22 of 57

22

Table Arrays

  • Elements in cell arrays are columns, each assigned a header name

  • same datatype per column (can be strings, numbers or other )
  • use table function to create such arrays
  • access columns using .

Example:

» T = table(var1, … , varN)

» T = table(Age,Height,Weight,'RowNames', LastName)

» disp (T.Age)

» meanHeight = mean (T.Height)

  • use readtable/writetable functions to read from/to the file

23 of 57

Next: Good Programming Practices

23

Questions?

24 of 57

find

24

  • find is a very important function
    • Returns indices of nonzero values
    • Can simplify code and help avoid loops

  • Basic syntax: index=find(cond)

» x = rand(1,100);

» inds = find (x>0.4 & x<0.6);

inds is a vector which contains the indices at which x has values between 0.4 and 0.6. This is what happens:

x>0.4 returns a vector with 1 where true and 0 where false x<0.6 returns a similar vector

& combines the two vectors using logical and operator find returns the indices of the 1's

25 of 57

Example: Avoiding Loops

25

length(x)

Loop time

Find time

100

0.01

0

10,000

0.1

0

100,000

0.22

0

1,000,000

1.5

0.04

  • Given x= sin(linspace(0,10*pi,100)), how many of the entries are positive?

Using a loop and if/else

count=0;

for n=1:length(x) if x(n)>0

count=count+1; end

end

Being more clever

count=length(find(x>0));

  • Avoid loops!
  • Built-in functions will make it faster to write and execute

26 of 57

Logical vectors

26

  • Logical vectors - consist of logical true/false values
  • For instance, to count how many vector elements are greater than 5:

» vect = [5 9 3 4 11];

» is_greater5 = vect>5;

» sum (is_greater5)

is_greater5 =

0 1 0 0 1

ans=

2

» vect = [5 9 3 4 11];

» is_greater5 = 0;

» for i = 1:5

» if vect(i)>5

» is_greater5 = is_greater5 + 1;

» end

» end

    • Longer and less efficient

27 of 57

Logical vectors

27

  • Logical vectors - consist of logical true/false values
  • Attention: 1s and 0s are logical, not integer!
  • To convert vector to logical use casting:

» a = logical ([1 0 1 0 1])

» b = int8(a)

a =

1×3 logical array

1 0 1 0 1

b =

1×3 int8 row vector

1 0 1 0 1

Use true-false instead or logical 0s to avoid confusion:

a = true false true false true

28 of 57

Logical vectors

28

  • Functions false and true create logical vectors and matrices

» a = false (2)

a =

0 0

0 0

» b = true (1,5)

b =

1 1 1 1 1

  • Logical function any(v) - returns true if at least one element of vector v is true
  • Logical function all(v) - returns true if all elements of vector v are true

» any ([1 0 1 0 1]) » all ([1 0 1 0 1])

ans = ans =

1 0

29 of 57

Vectorization

29

  • Avoiding loops
    • This is referred to as vectorization
  • Vectorized code is more efficient for MATLAB
  • Use indexing and matrix operations to avoid loops
  • For instance, to add 50 to each element of vector v:

v = randi([0 100],1,10);

for i=1:10

v(i) = v(i) + 50;

end

disp(v)

v = randi([0 100],1,10);

v = v + 50;

disp(v);

    • Efficient and clean

    • Slow and complicated

30 of 57

Vectorization

30

Sum up all elements of matrix A:

M = [1 , 2 , 1; 4 , 5, 2; 1, 3, 2];

sum_of_elements = 0;

number_of_rows = size (M , 1) ;

number_of_columns = size (M , 2) ;

% for each row

for row = 1: number_of_rows

% for each column

for column = 1: number_of_columns

sum_of_elements = sum_of_elements + M( row , column ) ;

end

end

disp ( sum_of_elements );

M = [1 , 2 , 1; 4 , 5, 2; 1, 3, 2];

sum_of_elements = sum(sum(M));

disp ( sum_of_elements );

    • Efficient and clean

    • Slow and complicated

31 of 57

Performance Measures

31

  • It can be useful to know how long your code takes to run
    • To predict how long a loop will take
    • To pinpoint inefficient code

  • You can time operations using tic/toc:

» tic % start timer

» My_code1;

» a= toc; % time 1

» My_code2;

» b = toc; % time 2

    • tic resets the timer
    • Each toc returns the current value in seconds
    • Can have multiple tocs per tic

32 of 57

Performance Measures

32

  • Sparse matrices

Sparse matrices occupy a lot of memory, leading to inefficient codes

Function sparse converts a full matrix into sparse form by squeezing out any zero elements. If a matrix contains many zeros, converting the matrix to sparse storage saves memory.

  • Example:

» A = zeros(10000); A(1,3) = 10; A(21,5) = pi;

» B = sparse(A);

» tic; A*A; toc % very slow

Elapsed time is 19.949076 seconds

» tic; B*B; toc % faster!

Elapsed time is 0.07464 seconds

33 of 57

Performance Measures

33

  • For more complicated programs, use the profiler

» profile on

    • Turns on the profiler. Follow this with function calls

» profile viewer

    • Displays gui with stats on how long each subfunction took

34 of 57

Next: Mini-Project

34

Questions?

35 of 57

35

HW5: Mini-project on Data Analysis

1. Download flight data for from the dataset alldelays.csv to access flight delays information provided by EUROCONTROL

2. Choose a research question from the list and clearly state it, add your group number in front of the corresponding questions

3. Extract the data you need for the corresponding airports/data periods from the dataset

4. Code the calculations to answer your research question

5. Illustrate the results using 2 different types of plots in Matlab

6. Send HW5 report to the examiner (Tatiana Polishchuk) in one email, containing:

  • Research question (from step 2)
  • Code covering the data extraction (step 3)
  • Code with calculations (step 4)
  • Code for two figures/plots created (step 5) and the figures in .jpeg or .png format attached

36 of 57

36

Advanced File Operations

EUROCONTROL Performance Review Unit (PRU): ATFM arrival flight delays

Data covers:

Years 2016-2024 (not much data for the pandemics years 2020-2021)

All European airports

37 of 57

37

Mini-project

EUROCONTROL Performance Review Unit (PRU): ATFM arrival flight delays

38 of 57

38

Advanced File Operations

Example: extract data subset from a data file

% clear Workspace

clear; clear all

% load flight delay data into the array variable delay_data

delay_data = readtable('alldelays.csv');

39 of 57

39

Advanced File Operations

Example: extract data subset corresponding for Munich airport

% extract delay data for the particular airport of Munich

flight_data = delay_data(strcmp(delay_data.APT_ICAO, 'EDDM'), :);

% save the delay data for Munich airport in a separate data file

writetable(flight_data, 'munich.csv');

40 of 57

40

Advanced File Operations

EUROCONTROL Performance Review Unit (PRU): ATFM arrival flight delays

Query examples:

  • What is the total number of flights delayed at Munich airport in July 2019?
  • What is the total amount of delays in minutes recorded at Munich airport in 2023?
  • For how many days in May 2023 the airport has experienced delays?
  • What were the main causes of Munich airport delays in 2023?
  • What is the share of the airport delays attributed to weather at Munich airport in 2023?
  • Compare the delay records (in total number of delayed flights and minutes of delays) at two airports of your choice, etc.

41 of 57

41

Advanced File Operations

Example: analyse the delays in Munich airport: code

% Clean the data to leave only the rows containing information about the delays

num_lines = size(flight_data,1)

count = 0;

for i=1:(num_lines-1)

if (flight_data.DLY_APT_ARR_1{i}>0)

count = count +1;

delay_data_munich(count,:) = flight_data(i,:);

end

end

writetable(delay_data_munich,'munich_delays.csv')

% now delay_data_clean contains only non-zero rows with delay information

% make calculations for the delays attributed to weather in Munich airport

42 of 57

42

Advanced File Operations

Example (cont.): analyse the delays in Munich airport code

% Calculate delays (in number of flights) Munich airport for all years

sum_delays = sum(delay_data_munich.FLT_APT_ARR_1)

% Calculate delays (in minutes) attributed to weather in Munich airport for all years

sum_weather_delays = sum(delay_data_munich.DLY_APT_ARR_W_1)

! Note: the data cells corresponding to minutes of delays is stored in “string” format in the delay datatable. This type is not appropriate for operation of summation!

There are 2 ways to resolve:

-use casting operation per cell in a loop applying str2double (it does not work for the whole column),

-or to make your life easy, open the datafile 'munich_delays.csv' (see previous slide) in Excel and format the corresponding column/cells to Number.

43 of 57

43

Plotting Figures

Example (cont.): analyse the delays in Munich airport code

%% Create a bar plot

figure;

bar([sum_delays, sum_weather_delays], 'FaceColor', [0.5 0.7 1]);

xticklabels({'Total Delays', 'Weather Delays'});

ylabel('Delay Minutes');

title('Munich Airport Delays');

grid on;

% Save figure

saveas(gcf, 'munich_delays_bar.png');

44 of 57

44

HW5: Mini-project on Data Analysis

More example figures

Fig1: Illustrate the proportions of different delay causes in airport of Munich in the year 2024

Fig2: Compare the impact of weather on delays at 3 airports for different years

Piechart: code Barplot: code

45 of 57

45

HW5: Mini-project on Data Analysis

*Tips: There are many plotting functions available in Matlab. Try different types, not only the ones I’ve showed in the example.

Very useful type of plotting statistics is boxplot (see an example below, but note that it does not correspond to the delay dataset).

Boxplot: code

46 of 57

Next: Mini-Project

46

47 of 57

47

HW5: Mini-project on Data Analysis

1. Download flight data for from the dataset alldelays.csv to access flight delays information provided by EUROCONTROL

2. Choose a research question from the list and clearly state it

3. Extract the data you need for the corresponding airports/data periods from the dataset

4. Code the calculations to answer your research question

5. Illustrate the results using 2 different types of plots in Matlab

6. Send HW5 report to the examiner (Tatiana Polishchuk) in one email, containing:

  • Research question (from step 2)
  • Code covering the data extraction (step 3)
  • Code with calculations (step 4)
  • Code for two figures/plots created (step 5) and the figures in .jpeg or .png format attached

48 of 57

48

Advanced File Operations

EUROCONTROL Performance Review Unit (PRU): ATFM arrival flight delays

Data covers:

Years 2016-2024 (not much data for the pandemics years 2020-2021)

All European airports

49 of 57

49

Mini-project

EUROCONTROL Performance Review Unit (PRU): ATFM arrival flight delays

50 of 57

50

Advanced File Operations

Example: extract data subset from a data file

% clear Workspace

clear; clear all

% load flight delay data into the array variable delay_data

delay_data = readtable('alldelays.csv');

51 of 57

51

Advanced File Operations

Example: extract data subset corresponding for Munich airport

% extract delay data for the particular airport of Munich

flight_data = delay_data(strcmp(delay_data.APT_ICAO, 'EDDM'), :);

% save the delay data for Munich airport in a separate data file

writetable(flight_data, 'munich.csv');

52 of 57

52

Advanced File Operations

EUROCONTROL Performance Review Unit (PRU): ATFM arrival flight delays

Query examples: full list

  • What is the total number of flights delayed at Munich airport in the summer 2019?
  • What is the total amount of delays in minutes recorded at Munich airport in 2023?
  • For how many days in May 2023 the airport has experienced delays?
  • What were the main causes of Munich airport delays in 2024?
  • What is the share of the airport delays attributed to weather at Munich airport in 2023?
  • Compare the delay records (in total number of delayed flights and minutes of delays) at two airports of your choice, etc.

53 of 57

53

Advanced File Operations

Example: analyse the delays in Munich airport: code

% Clean the data to leave only the rows containing information about the delays

num_lines = size(flight_data,1)

count = 0;

for i=1:(num_lines-1)

if (flight_data.DLY_APT_ARR_1{i}>0)

count = count +1;

delay_data_clean(count,:) = flight_data(i,:);

end

end

writetable(delay_data_clean,'munich_delays.csv')

% now delay_data_clean contains only non-zero rows with delay information

% make calculations for the delays attributed to weather in Munich airport

54 of 57

54

Advanced File Operations

Example (cont.): analyse the delays in Munich airport code

% Calculate delays (in number of flights) Munich airport for all years

sum_delays = sum(delay_data_clean.FLT_ARR_1)

% Calculate delays (in minutes) attributed to weather in Munich airport for all years

delay_data_clean_double = readtable('munich_clean.csv');

sum_weather_delays = sum(delay_data_clean_double.DLY_APT_ARR_W_1)

! Note: the data cells corresponding to minutes of delays is stored in “string” format in the delay datatable. This type is not appropriate for operation of summation!

There are 2 ways to resolve:

-use casting operation per cell in a loop (it does not work for the whole column),

-or to make your life easy, open the datafile 'munich_delays.csv' (see previous slide) in Excel and format the corresponding column/cells to Number.

55 of 57

55

HW5: Mini-project on Data Analysis

EXAMPLE:

Fig1: Illustrate the proportions of different delay causes in airport of Munich in the year 2024

Fig2: Compare the impact of weather on delays at 3 airports for different years

Piechart: code Barplot: code

56 of 57

56

HW5: Mini-project on Data Analysis

*Tips: There are many plotting functions available in Matlab. Try different types, not only the ones I’ve showed in the example.

Very useful type of plotting statistics is boxplot (see an example below, but note that it does not correspond to the delay dataset).

Boxplot: code

57 of 57

THANK YOU!

57

Questions?