Using finite differences to find the displacement (x,y) of a projectile including the effect of drag

When i enter the code below I get the following error:

Index exceeds the number of array elements (2).
Error in EGB211_CLA (line 37)
    XxNow = XxN(i-1) ;

Would anyone know how to fix this issue? My code is below as well

close all
clc
 
%% Problem Parameters 
 
m = 0.05 ;  % mass (kg) = 50x10^-3
angle = 25 ; % initial angle (degrees)
th0 = angle*(pi/180) ; % initial angle (radians)
Xx0 = 0 ; % initial horizontal displacement (m)
Xy0 = 50 ; % initial vertical displacement (m)
V0 = 50 ; % initial velocity magnitude (m/s)
Vx0 = 41.315 ; % initial horizontal velocity (m/s)
Vy0 = 21.1309 ; % initial vertical velocity (m/s)
r = 0.03 ; % radius (m)
A = pi*(r)^2 ; % cross sectional area of projectile (m^2)
Cdrag = 0.80 ; % drag coefficient (unitless)
rho = 1.183 ; % air density (kg/m^3)
g = 9.81 ; % gravity (m/s^2)
 
 
%% Numerical Solution (With Drag) 
 
t = 10 ;
dt = 0.01 ; % Time step (seconds)
Niter = floor(10/0.01) ; % total number of intervals (t/dt)
 
XxN(1) = Xx0 ; 
XxN(2) = Xx0+(dt*Vx0) ; 
    
XyN(1) = Xy0 ; 
XyN(2) = Xy0+(dt*Vy0) ; 
 
for i = 3:Niter
    
    XxNow = XxN(i-1) ; 
    XxPrev = XxN(i-2) ;
    
    XyNow = XyN(i-1) ; 
    XyPrev = XyN(i-2) ;
    
    fx = @(XxNext)m*((XxNext-(2*XxNow)+XxPrev)/(dt^2)+(0.5*Cdrag*rho*A)...
        *(sqrt(((XxNow-XxPrev/dt)^2)+((XyNow-XyPrev/dt)^2))))*(XxNext-XxPrev/2*dt) ;
    
    fy = @(XyNext)m*((XyNext-(2*XyNow)+XyPrev)/(dt^2)+(0.5*Cdrag*rho*A)...
        *(sqrt(((XxNow-XxPrev/dt)^2)+((XyNow-XxPrev/dt)^2))))*(XyNext-XyPrev/2*dt)+(m*g) ;
    
    XxNext = fzero(fx,XxNow) ; 
    xxN(i) = XxNext ; 
    
    XyNext = fzero(fy,XyNow) ; 
    xyN(i) = XyNext ; 
    
end 

Answer 1

[KSSV]

There was a type error in updating XxN.....it was typed as xxN instead of XxN.....find the below corrected code:

close all
clc

%% Problem Parameters 

m = 0.05 ;  % mass (kg) = 50x10^-3
angle = 25 ; % initial angle (degrees)
th0 = angle*(pi/180) ; % initial angle (radians)
Xx0 = 0 ; % initial horizontal displacement (m)
Xy0 = 50 ; % initial vertical displacement (m)
V0 = 50 ; % initial velocity magnitude (m/s)
Vx0 = 41.315 ; % initial horizontal velocity (m/s)
Vy0 = 21.1309 ; % initial vertical velocity (m/s)
r = 0.03 ; % radius (m)
A = pi*(r)^2 ; % cross sectional area of projectile (m^2)
Cdrag = 0.80 ; % drag coefficient (unitless)
rho = 1.183 ; % air density (kg/m^3)
g = 9.81 ; % gravity (m/s^2)


%% Numerical Solution (With Drag) 

t = 10 ;
dt = 0.01 ; % Time step (seconds)
Niter = floor(10/0.01) ; % total number of intervals (t/dt)

XxN = zeros(Niter,1) ; 
XyN = zeros(Niter,1) ; 

XxN(1) = Xx0 ; 
XxN(2) = Xx0+(dt*Vx0) ; 

XyN(1) = Xy0 ; 
XyN(2) = Xy0+(dt*Vy0) ; 

for i = 3:Niter
  
  XxNow = XxN(i-1) ; 
  XxPrev = XxN(i-2) ;
  
  XyNow = XyN(i-1) ; 
  XyPrev = XyN(i-2) ;
  
  fx = @(XxNext)m*((XxNext-(2*XxNow)+XxPrev)/(dt^2)+(0.5*Cdrag*rho*A)...
  *(sqrt(((XxNow-XxPrev/dt)^2)+((XyNow-XyPrev/dt)^2))))*(XxNext-XxPrev/2*dt) ;
  
  fy = @(XyNext)m*((XyNext-(2*XyNow)+XyPrev)/(dt^2)+(0.5*Cdrag*rho*A)...
  *(sqrt(((XxNow-XxPrev/dt)^2)+((XyNow-XxPrev/dt)^2))))*(XyNext-XyPrev/2*dt)+(m*g) ;
  
  XxNext = fzero(fx,XxNow) ; 
  XxN(i) = XxNext ; 
  
  XyNext = fzero(fy,XyNow) ; 
  XyN(i) = XyNext ; 
  
end 

[Finnien Mitchell]

Ok that works but now matlab will only continue to run and won't stop unless I pause it and tell it to quit debugging which is an option I have not seen before. Would you know how to fix this?

[KSSV]

I suspect you made a breaking point ...remove that point.

[Finnien Mitchell]

A breaking point? I entered the same code as above?

[KSSV]

Code should run fine....you close your file and open once....you are running code in debug mode I think.

[Finnien Mitchell]

i've closed the file and opened it and its still not working. Maybe I have a breakpoint in my code somewhere??

[KSSV]

The code is working fine for me......Breakingpoint can be easily known...it will appear as a red dot on the left of the editor.

[Finnien Mitchell]

oh ok thankyou