Newton Method for solving 2nd order nonlinear differential equation

Hello!

I have the following code for solving a 2nd order nonlinear differential equation, which a mass-spring-damper system with an external forcing f.cos(omega*t) and the nonlinear term in alpha*x^3.

When I attempt to run the code it gives me the following error:

Error using vertcat

Dimensions of arrays being concatenated are not consistent.

Error in NewtonTry1 (line 26)

F = [F1(y) ; F2(y)];

My script is the following:

clc 
clear all
 
% defining parameters and variables
f = 1;
n = 32;
t = 0:(n-1);
m = 1;
k = 1;
c = 0.01;
alpha = 0.1; 
omega = 2*pi/n;

% initial guess 
y = [0 ; 0];

% defining 2nd order ODE as a system of 2 1st order ODEs
F1 = @(y) y(2);
F2 = @(y) (-c/m)*y(2)-(k/m)*y(1)-(alpha/m)*y(1)^3+(f/m)*cos(omega*t);
%dydt = dydt[:];
% F1 = dydt(1);
% F2 = dydt(2);

for i = 1:10
    %F = zeros(size(y));
    F = [F1(y) ; F2(y)];
    
    J11 = deriv(F1,y,y(1),1);
    J12 = deriv(F1,y,y(2),2);
    J21 = deriv(F2,y,y(1),1);
    J22 = deriv(F2,y,y(2),2);
    
    Jac = [J11 J12;
           J21 J22];
       
       y = y - inv(Jac)*F;
    
end

And the derivative user-defined function written as such:

% creating a user-defined function for calculation of derivatives
function Jac = deriv(F,y,yi,i)
% yi = which element of y
% i = element number

% numerical differentiation f'(x) = (f(x+e)-f(x))/e
e = 0.0001;
y1 = y;
y2 = y;
y1(i) = yi;
y2(i) = yi+e;
Jac = (F(y1)-F(y2))/e;

end
  

How do I go about fxing this error? Any suggestions and advice would be very much appreciated!

Answer 1

[Divyam]

Hi @Vedika Pandya,

You are concatenating two arrays of different sizes and hence are facing this issue, since t is an array, F2(y) is a 1x32 matrix while F1(y) is simply an integer 0. Try taking individual values of t to obtain a value of F2(y) which can be concatenated with F1(y).