Write an algorithm describing how this would be implemented to solve for given h, , and as inputs.

Try solving the following ODEs using ode45.

1.

% Insert code here

2.

% Insert code here

1

3.

% Insert code here

Solving Systems of ODEs using ode45

If the system your are solving for is a higher order ODE, or can’t be fully described as a single ODE, then we must solve by treating it as a system of ODEs. Solving a system of ODEs using ode45 is the same as what we saw previously, except the inputs are now vectors.

The first input, , is now a column vector containing equations for y’. We now write it in the

form . The input vairable y is a column vector as well, and so to access the individual dependent variables, we call their location in the y vector: y(1) and y(2).

The second input, , is unchanged form what we saw previously.

The third input, , is now a column vector containing the initial values for all the variables.

is still a vector containing t values, but is now a matrix containing y values, where each column is a different solution vector.

Example:

ODEFUN = @(t,y) [y(2)-y(1);2*y(1)-y(2)]; TSPAN = [0 5]; Y0 = [1;0.5]; [TOUT,YOUT] = ode45(ODEFUN,TSPAN,Y0); plot(TOUT,YOUT) legend(‘y_1′,’y_2’)

Now try the following examples.

1.

% Insert code here

2.

% Insert code here

3.

% Insert code here

2

4.

% Insert code here

Writing Your Own ODE Solver

When MATLAB is readily available, using will be the simplest and most accurate way to numerically solve ODEs. However, if MATLAB is unavailable, it becomes important to understand how to write your own ODE solver. This will also give you some more practice with writing algorithms. In last week’s lecture, we went through three common numerical methods for solving ODEs, which were Euler’s method, modified Euler’s method, and Runge-Kutta (RK4). We saw how to code the Euler method and modified Euler method in class and now you will implement the RK4 method.

RK4 Method: the RK4 method, just like Euler’s method, solves first order ODEs in the form . The method solves the ODE by evaluating the following equations.

Where the k values are defined through the following equations:

We can see that the RK4 method will require four inputs. The first three are the same as ode45, in that we need to know the ODE we are solving, the time span and the initial condition. The fourth input required will be the step size, h. We would like for the function to create two outputs, TOUT and YOUT. The process of applying RK4 is iterative, meaning we should be able to describe it in MATLAB using a for loop or a while loop. Now try and write the following numerical solvers.

1. Write an algorithm for the RK4 method. 2. Convert the algorithm to code that computes the RK4 method. You can submit your code on AMS

(practice activities) to check that it’s correct. 3.

Solve the ODE using your RK4 function, given , , y(0) = 1.

4. Compare your solution to the Euler, modified Euler, and ode45 solutions (these are all available to download from the week 8 lecture files (example1.m, example2.m, example3.m).

5. *Bonus (time permitting): experiment with the accuracy of each method as the step size is changed. Try seeing if you can do this using a for loop.

3

Portfolio 3 Part C – competing species

In Central Queensland, cattle farmer Sam is trying to grow sufficient grass for his cattle to feed on. However, as his grass grows, the kangaroos start to gather and eat Sam’s grass! To combat this, Sam plants additional grass on his land that, in turn, attracts even more kangaroos. To understand the dynamic of these competing species (the grass (g(t)) and kangaroos (k(t)), we can use the following system of ODEs.

where and are growth constants.

1) Write out explicitly the matrix-vector representation of this system.

Clearly define the elements

2) Write an algorithm describing how this would be implemented to solve

for given h, , and as inputs.

3) Implement this algorithm in MATLAB to solve for given and you can assume that initially there are no Kangaroos and 1083ha of consumable grass available.

4) Plot the solutions (in an appropriate manner). Explain why this model is or isn’t appropriate (is there a point where it becomes unrealistic?)

4

3parts-portfolio.docx

notes-partA.pdf

Algorithms and control structures: part a

This week we will introduce algorithm development and two control structures: if statements and for loops. Next week, you will be introduced to while loops and nested conditional statements. The portfolio tasks will be combined into a single week.

Algorithms

One of the hardest transitions to make when first learning to program is how to ‘think like a computer’. Computers can do simple tasks extremely fast, but must, at some point, be given specific instructions on what to do.

A computer performs all of its tasks through the application of an algorithm. An algorithm is a set of instructions that can be performed one by one. They are as basic as possible (although can add together to produce complex behaviour). Each step in an algorithm should be irreducible enough to remove any ‘open’ interpretation – that is, given the same inputs, the same result should be obtained.

In general, when writing an algorithm, it is a good idea to provide four key parts: Input, Setup, Process, and Output. In additional to this, the algorithm needs to terminate.

Input: should contain all values given to the algorithm by a user, except where prompts for values are included in the algorithm process.

Setup: involves arranging the initial state of the system, such as initialising a sum variable to zero.

Process: is the essential part of the algorithm, describing how the output is obtained.

Output: should list things that specifically are given by the algorithm, the set of results that you should obtain when following the algorithm.

Terminate: Finally, all Algorithms must come to an end.

Activity 1

Imagine that you have access to a “human computer”. This HC is capable of taking instructions in English (so long as the instructions are simple ones that are clear and concise), and can perform basic human actions, such as “put on shirt”. The HC understands concepts such as closeness, temperature, time, and the various weather conditions.

1. Working alone, write an algorithm for the process of getting dressed in the morning.

2. Compare your algorithm to other members in your group. How different are they? Is there any room for interpretation?

From Algorithm to Program

There are three forms that we could write instructions in.

Algorithm: is a series of instructions in human language, providing a clear and concise way for a human to read and understand the process involved. Algorithms are often written to suit the level of

the human that is expected to read it, and thus may include complex instructions that it is assumed the reader knows, such as “obtain the roots of the quadratic”.

Pseudocode: is a set of instructions written in a generic language. That is, it is written to still be easily followed by a human, but is formatted like most programming languages. It often uses simple human instructions such as “Set x=2”, “Repeat Until x>5”, or “Stop When err<0.001”. When writing pseudocode, it is often beneficial to format your pseudocode with mathematical notation, such as e^x, rather than computer notation, such as exp(x).

Code: is written for a specific computer language – in this course, you will write your code for MATLAB. We will write our code in “m-files”, which should be saved to the current working directory. Code in MATLAB comes in two forms – scripts and functions. Scripts are simply a series of instructions that are fed into MATLAB as though they were typed into the command prompt. Functions take inputs and are called in the same way as built-in functions.

Programming: If Statements

An ‘if’ statement is like a branch in the process. When the computer hits the branch, it goes in one of two directions, depending on whether a condition is true or false. It does this by first testing the condition using ‘if’, before giving one option. Multiple branches can be expressed by using ‘elseif’ and ‘else’.

For the examples below, run for x values of 1,5 and 10 and see what happens.

Example 1: If statement

x = 1; if x > 5 x = x + 2 end

Example 2: Elseif statement

x = 1; if x == 10 x = 0; elseif x > 4 x = x + 2 end

Example 3: Else statement

x = 1; if x >= 5 x = x + 2 else error(‘number is less than 5’) end

Activities

1. Write an if statement that dispays the value of a variable (x) if the value is less than or equal to zero.

% Write your code here

2. Write an algorithm that finds the larger of two numbers (a,b). Consider positives to always be larger than negatives (Eg. 3 is larger than -5, -1 is larger than -4). Convert this to a MATLAB function called max2.

% Write your code here

3. Write an algorithm that finds the largest of the four numbers (a, b, c, d). Convert this to a MATLAB function called max4. hint: you do not need to check every single combination

% Write your code here

4*. Alter the previous code so that it counts the number of times the maximum number appears.

% Write your code here

Programming: For Loops

A ‘for’ loop repeats a process a set number of times, but changes one variable (called the counter) each time. They are particularly useful when dealing with any repepitive process. The range of values that the counter takes on is given at the start of the for loop.

Run the examples below to see how for loops work.

Example 1: Basic For Loop

The following for loop displays the whole numbers from 1 to 10.

for i = 1:10 disp(i) end

Example 2: Triangular Numbers

The following for loop determines the Nth triangular number.

N = 10; %input triangular = 0; %setup – initial value required for i = 1:N triangular = triangular + i end disp(triangular) %output

Example 3: Triangular Numbers with Storage

The following for loop stores all the triangular numbers up to the Nth number within a vector.

N = 10; %input triangular(1) = 1; %setup – initial value required for i = 2:N triangular(i) = triangular(i-1) + i end disp(triangular) %output

Activities

1. Write a code using a for loop that displays the numbers 1 to 10.

% Write your code here

2. Write a code that sums the squares of the numbers 1 to 100.

% Write your code here

3. Write a code stores the sum of squared numbers 1 to 100.

% Write your code here

4. Write a code that stores the squared even numbers from 2 to 20.

% Write your code here

Numerical differentiation

In class you saw some code that found the difference between cells in a vector: we can build on this to numerically differentiate a function. We will utilise a central difference approximation

where h is the size of the increment of the independent variable (this is oftern refered to as the step size).

Write an algorithm that calculates the numerical derivative of a function. Convert this algorithm to code, the function has been started below for you. Make sure you plot both the function and it’s derivative to make sure you have the right answer.

f = @(x) x.^2+2*x x = 0:0.01:3; y = f(x); % Write your code here

Portfolio Task 3 Part a.

Submission instructions:

1. This portfolio should be submitted as a single pdf file titled n#######Portfolio3.pdf 2. The assignment is DUE at 11:59pm Friday 11th May and should be submitted through Blackboard.

The IT helpdesk closes at 10pm, and so if you leave your submission to the last minute and something happens then you may not be able to recieved any assistance. Late submissions will recieve a mark of zero.

Algorithms are a great way to visualise a mathematical method without worrying about the algebra. We’d highly recommend writing algorithms (and pseudocode) to help you with your study for the exam. The purpose of this assignment is to provide you an opportunity to try writing your own algorithm on a previous method learned in this class.

1) Write an algorithm that solves a second order non-homogenous ODE (with boundary conditions) using the method of undeterimined coefficients. (Note, this algorithm should include between 4 and 16 steps).

2) From your algorithm in question 1), write a pseudocode that describes how to solve the method of undertermined coefficients.