ODE45 not updating a variable of nonlinear 2nd order diff.eq system

Hello everyone, I was trying to solve a simple nonlinear system of diff.eqn (piece-wise linear) that involves coulomb friction. The problem I face is x(3) does not update after passing if condition... Although when i debug it, I can see that its varying accordingly... but at the end get a vector filled with a constant number(initial condition x3(0).) Could you please guide me how to cope with this issue???

function [dxdt]=func6(t,x)

C=10;
k=5;
kd=10;
F=25;
mu=0.5;
Nf=15;
m=1;
wn=sqrt(k/m);
w=wn;
fE=F*sin(w*t);
dxdt=zeros(size(x));
x1=x(1);
x2=x(2);
x3=x(3);
dxdt(1)=x(2);

%

nonlinear force calculation function

F_NL=kd*(x(1)-x(3));

if F_NL>=mu*Nf

     F_NL=mu*Nf;

   x(3)=x(1)-(mu*Nf)/kd;
  end

if  F_NL<=-mu*Nf

      F_NL=-mu*Nf;

      x(3)=x(1)+(mu*Nf)/kd;
  end

dxdt(2)=1/m*(fE-F_NL-C*x(2)-k*x(1));

end

%%%and am using this script to run my ode calculations...

h=0.01;

tspan = 0:h:10;

initialvalues=[0 0 0];

[t,x]=ode23(@func6,tspan,initialvalues);

Answer 1

[Aquatris]

The reason is simple. Your dxdt(3) is always equal to 0. You do not have an equation for dxdt(3) in your function other than the first zeros(size(x)) assignment. Since the initial value is also zero, x(3) stays 0 at all times.

By any chance, those if loops supposed to be for dxdt(3) not x(3)?

[Res]

Thank you for replying, yes I realized so , but with this system i don,t have any eqn for dxdt(3) .. My point is , can I save values of x(3) at each iteration?

[Aquatris]

It is usually not recommended but one way is to use global variable. Instead of using x(3), create a global variable x3 and assign x3(length(x3)+1) whenever you need to store a new value to x3.

[Res]

Could you please elaborate little bit, am quite beginner on matlab. As my x(3) is not updating in the initial conditions (because i don,t have eqn in the form of dxdt(3)...), so at each iteration F_NL=kd*(x(1)-x(3)); runs with F_NL=kd*(x(1)-x(0)); which is wrong... how can I save this value x(3)=x(1)+(mu*Nf)/kd and recall it in F_NL=kd*(x(1)-x(3)); Thanks

[Aquatris]

Something like this, it needs some modifications though;

The function;

function [dxdt]=func6(t,x)
 % specify the global variables
 %   the function needs to use
 global x3 counter t_x3 

   C=10;
   k=5;
   kd=10;
   F=25;
   mu=0.5;
   Nf=15;
   m=1;
   wn=sqrt(k/m);
   w=wn;
   fE=F*sin(w*t);
   dxdt=zeros(size(x));
   x1=x(1);
   x2=x(2);

   dxdt(1)=x(2);

nonlinear force calculation function

   F_NL=kd*(x(1)-x3(counter));
   counter = counter +1;% counter to store x3 values
   t_x3(counter) = t;   % time for variable x3, ode is weird sometimes
                        %    and the size of the time output of ODE did 
                        %    not match the  x3 variable
   if F_NL>=mu*Nf
       F_NL=mu*Nf;
       x3(counter)=x(1)-(mu*Nf)/kd;
   elseif  F_NL<=-mu*Nf
       F_NL=-mu*Nf;
       x3(counter)= x(1)+(mu*Nf)/kd;
   else
       x3(counter) = 0; % I assigned 0 but what happens here 
                        %   is up to your problem             
   end
   dxdt(2)=1/m*(fE-F_NL-C*x(2)-k*x(1));
  end

And the script;

clear,clc
% specify the global variables
global x3 counter t_x3
x3 = 0;
counter = 1;

h=0.01;
tspan = 0:h:10;
initialvalues=[0 0 ]';
[t,x]=ode45(@func6,tspan,initialvalues);

plot(t,x,t_x3,x3,'.')

Answer 2

[Res]

Thank you, really appreciated.