# Where/How exactly do I create an ODE Event Location Function?

51 views (last 30 days)

Show older comments

Hello all,

Here's the thing:

I have a main file (MainSimulation.m), which basically creates a menu in which I can choose a type of marine Manoeuvre to perform, given variables in other two files (data.m and wind.m), which are the vessel-related parameters and the wind conditions + wind coefficient (forces/moment) calculations.

Now, when I choose a Manoeuvre, I must run ode45 to be able to calculate the overall accelerations (F=m.a), like so:

option = menu('Choose a Manoeuvre:','Straight Run','Zig-Zag Manoeuvre','Exit Menu');

%% Straight Run

if option == 1 % Choosing Straight Run

Vs = 10 * (1852 / 3600); % [m/s]-Ship approach speed (reasonable value)

a_1 = 300;

t = (0:a_1); % Manouevre time

y0 = [Vs 0 0 0 0 0 0]; % Initial condition - surge valocity is equal to the approach speed

[T,Y] = ode45 ('StraightRun', t, y0); % Calls function

drift_angle = -atan(Y(:,2) ./ Y(:,1));

In which Y will be a 7 collumn double array [Y(:,1) = speed in X; Y(:,2) = speed in Y; the other collumns are irrelevant rn]

Well the thing is, given the fact that this is a wind-powered vessel (obtained from data.m) I have come to the conclusion that there is a problem with the ODE function (any of them) when the speed (Vs) of the vessel comes to 0, because there are several parameters (regarding the definition of forces/moments) that directly depend on the vessel's speed, and the following Warning appears, if the ODE45 is not done running, given t (which I assume to be due to the vessel coming to a full stop and there still being some integration steps left):

Warning: Failure at t=3.091528e+02. Unable to meet integration tolerances without reducing

the step size below the smallest value allowed (9.094947e-13) at time t.

> In ode45 (line 351)

Therefore, I wish to run the ODE45 function, as stated, until Vs = 0.

I have seen that there is a thing called ODE Event Location which allows for, if a certain variable/function comes to 0, the ODE will stop integrating, and return the values already calculated.

However, I'm not really sure WHERE to insert the following function [and I'm a little confused as to what "isterminal" and "direction" are]:

function [position,isterminal,direction] = EventsFunction(t,y)

position = sqrt(y(1)^2+y(2)^2); % The value that we want to be zero

isterminal = 0; % Halt integration

direction = 0; % The zero can be approached from either direction

end

Should I insert the above function (EventsFunction) in the same file as the menu is, prior to it, obviously, or should I insert it inside the 'StraightRun' file?

Because I have tried inserting the function "EventsFunction" just before calling ode45 (after the definition of y0), but MATLAB gives me the following error:

Error: File: MainSimulation.m Line: 32 Column: 1

Function definition are not supported in this context. Functions can only be created as local

or nested functions in code files.

And if I insert the function BEFORE/outside the option = menu (....), the following error appears:

Error: File: MainSimulation.m Line: 29 Column: 1

Function definitions in a script must appear at the end of the file.

Move all statements after the "EventsFunction" function definition to before the first local

function definition.

ANY Help provided will be MUCH appreciated!

Nuno Nunes

##### 0 Comments

### Accepted Answer

Torsten
on 13 Aug 2022

if option == 1 % Choosing Straight Run

Vs = 10 * (1852 / 3600); % [m/s]-Ship approach speed (reasonable value)

a_1 = 300;

t = (0:a_1); % Manouevre time

y0 = [Vs 0 0 0 0 0 0]; % Initial condition - surge valocity is equal to the approach speed

opt = odeset('Events',@EventsFunction);

[T,Y] = ode45 ('StraightRun', t, y0, opt); % Calls function

drift_angle = -atan(Y(:,2) ./ Y(:,1));

and place the function where MATLAB can access it.

function [position,isterminal,direction] = EventsFunction(t,y)

position = y(1)^2+y(2)^2; % The value that we want to be zero

isterminal = 1; % Halt integration

direction = 0; % The zero can be approached from either direction

end

##### 15 Comments

Jan
on 31 Aug 2022

@Nuno Nunes: What does "ode45 is giving discontinuous values" mean?

If the value of an event function is below 0 in one step and above 0 in the next step, ODE45 interpolates the step to find the time, where 0 is hit exactly.

"from my experience, it does not." - please post some code to explain this. Of course the show event function is valid and should be detected reliably. Maybe reducing the maximum step size helps, if the event function is > 0 at both steps, but was < 0 between the time points.

### More Answers (2)

Walter Roberson
on 14 Aug 2022

You are using

[T,Y] = ode45 ('StraightRun', t, y0, opt); % Calls function

specifying the function to be executed as a character vector. We recommend you stop doing that and start using function handles. Besides being more efficient, function handles are able to refer to local functions defined in the same file.

So you can create a single file

option = menu('Choose a Manoeuvre:','Straight Run','Zig-Zag Manoeuvre','Exit Menu');

%% Straight Run

if option == 1 % Choosing Straight Run

Vs = 10 * (1852 / 3600); % [m/s]-Ship approach speed (reasonable value)

a_1 = 300;

t = 0:a_1; % Manouevre time

y0 = [Vs 0 0 0 0 0 0]; % Initial condition - surge valocity is equal to the approach speed

opt = odeset('Events',@EventsFunction);

[T,Y] = ode45 (@StraightRun, t, y0, opt); % Calls function

drift_angle = -atan(Y(:,2) ./ Y(:,1));

end

function [position,isterminal,direction] = EventsFunction(t,y)

position = y(1)^2+y(2)^2; % The value that we want to be zero

isterminal = 1; % Halt integration

direction = 0; % The zero can be approached from either direction

end

function dy = StraightRun(t, y)

whatever appropriate

end

I have no idea what your data.m is about. You are not defining anything that needs to be passed into StraightRun.

If StraightRun needs additional data, see http://www.mathworks.com/help/matlab/math/parameterizing-functions.html

##### 2 Comments

Walter Roberson
on 15 Aug 2022

Jan
on 15 Aug 2022

As Walter has written already, use:

[T,Y] = ode45 (@StraightRun, t, y0);

% Instead of the old style:

% [T,Y] = ode45 ('StraightRun', t, y0);

Providing the function as CHAR vector is deprecated for 20 years now and "modern" function handles are the better choice.

Then your problem does not concern only the event function, but all functions. The rules are:

- A function can be stored in a specific m-file, which ist stored in a folder inside Matlab's path. See doc addpath and doc pathtool.
- A function can by stored as local subfunction inside an m file. If this file is a "script", the functions must be appended at the bottom. Scripts are m-files, which do not start with the keyword "function".
- Nested functions are stored inside another function. They share variables with the enclosing function.

Examples:

% m-file: myScript.m

disp('This is a script');

a = myFunc(2);

% Local subfunction:

function out = myFunc(in)

out = in * 3;

end

% m-file: myFunction.m

function myFunction % A function without inputs

disp('This is a script');

a = myFunc(2);

end % <- required to close the function, if any other function uses "end" also

% Local subfunction - only visible inside this file:

function out = myFunc(in)

out = in * 3;

end

% m-file: myFunction.m

function myFunction % A function without inputs

disp('This is a script');

a = myFunc(2);

% Nested function - only visible in enclosing function

function out = myFunc(in)

out = in * 3;

end

end % <- required to close the function, if any other function uses "end" also

% m-file: myFunction.m

function myFunction % A function without inputs

disp('This is a script');

a = myFunc(2);

end

% Another m-file: myFunc.m - visible to all other functions and scripts

function out = myFunc(in)

out = in * 3;

end

##### 0 Comments

### See Also

### Categories

### Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!