Using very time dependent functions with ode suite in MATLAB

Hello everyone,

I am trying to solve an ordinary differential equation that has a function in it that only activates when T=0.05 (essentially a dirac impulse function). I need o figure out a way to ensure that even with an adaptive time step this function will trigger. Below I've attached my ode function and my script to run the file. Any help would be greatly appreciated.

(Please note that I do not use the actual dirac function in the ode because it causes values to go to NaN because of the linear algebra involved. My ode function:

function [Xddot]=mdof(T,X)

%System Parameters
m=30; %kg
M=200;  %kg
a=3;  %m
b=1;  %m
k1=400*1000;  %N/m
k2=k1;
k3=100*1000;  %N/m
k4=k3;
I=200;  %kg*m^2
C1=3200;  %Ns/m
C2=C1;
C3=800;  %Ns/m
C4=C3;

M2=[(M*b^2+I)/(a+b)^2 (M*a*b-I)/(a+b)^2      0    0;...
   (M*a*b-I)/(a+b)^2 (M*a^2+I)/(a+b)^2      0    0;...
            0                 0             m    0;...
            0                 0             0    m];

C=[  C1   0   -C1   0 ;...
     0    C2   0   -C2;...
    -C1   0  C1+C3  0 ;...
     0   -C2   0  C2+C4];

K=[  k1    0    -k1     0 ;...
     0    k2     0     -k2;...
    -k1    0  (k1+k3)   0 ;...
     0    -k2    0   (k2+k4)];

   Xdot=[X(5);X(6);X(7);X(8)];
   XX=[X(1);X(2);X(3);X(4)];
   %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
   F=[0;0;(T==0)*1e100;((T-0.05)==0)*1e100]
   %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
   T
   Z=(F-C*Xdot-K*XX);

   Xddot2=M2\Z;

   Xddot=[Xdot;Xddot2];
  end

Here the part set off is where I am having the problem. The adaptive time step jumps over T=0.05 and thus never triggers the input. If i try to use a error tolerance though, such as

(abs(T-0.05)<=1e-2)

I either end up with the same problem of the function never triggering or a problem of the function triggering numerous times because the time step slows way down to handle the high input.

Here's the call function for if you want to run it.

clc, clear, close all
%%Numerical Method
IC=zeros(1,8);

t=0:0.0001:5;
options=odeset('RelTol',1e-6,'AbsTol',1e-6);
[T,X]=ode45('mdof',t,IC,options);

figure
plot(T,X(:,1),'r-')
hold on
plot(T,X(:,2),'k-')
plot(T,X(:,3),'g-')
plot(T,X(:,4),'-')
title('\textbf{Numeric Analysis Results}','interpreter','LaTeX','FontSize',16)
xlabel('Time, s','interpreter','LaTeX','FontSize',14)
ylabel('Displacement, m','interpreter','LaTeX','FontSize',14)
legend('x_1','x_2','x_3','x_4')
hold off

Thanks so much for your help!!

Answer 1

[Jan]

Matlab's ODE integrators cannot handle discontinuities, see http://www.mathworks.com/matlabcentral/answers/59582#answer_72047 .

The only stable and reliable method is to integrate until the jump, apply the changes and proceed with a new integration.