I have the following nonlinear system and accompanying ODE:


My solution in MATLAB:

% Solving with ODE45
dt = 10^-1; % [s]
T = 100; % [s]
N = 2; % number of masses
X_0 = zeros(2*N,1);
t_span = [0:dt:T];
m_i = 1; % [kg]
k_i = 2; % [N/m]
c_i = 1; % [kg/s]
gamma_i = 1; % [N/m^3]
f_i = 1; % [N]
f = @(t) f_i.*sin(alpha.*t).*cos(beta.*t); % [N]
m = m_i.*ones(N,1); % [kg]
k = k_i.*ones(N+1,1); % [N/m]
c = c_i.*ones(N+1,1); % [kg/s] 
gamma = gamma_i.*ones(N+1,1); % [N/m^3]
alpha = 0.2;
beta = 0.15;
K = karray(k,gamma,N);
M = marray(m);
C = carray(c,N);
F = @(t) farray(f,t,N);
options = odeset('Mass',M,'RelTol',1e-3,'AbsTol',1e-3);
fun = @(t,X)odefun(X,K,C,M,F(t),N); % Defines the ODE function and it's input matrices
[t_ode,X_answer] = ode15s(fun,tspan,X_0,options); % Calls the ODE function (see end of code), sets the time and axial conditions. Outputs time and displacements.

The functions karray(),carray(),marray(),farray() return array representations of the system parameters while odefun() returns an array of the EOM’s displacements and velocities. I had set the initial displacement and velocity of each mass to zero as my initial values.

I get the error “Error using matlab.internal.decomposition.builtin.luSolve.” when trying to solve with ode15s() or ode23t(). Is this a memory issue or am I not calling these functions correctly?