ode45 invalid parameter 981/100

clear;
animate = 1;    % 1 = animate, 0 = no animation

syms m L g y_0 omega t y_ddot theta theta_ddot
y = y_0*sin(omega*t)

y = 

eqn(1) = y_0*sin(omega*t) -m*g == m*y_ddot;
eqn(2) = 0.5*L*y_0*sin(theta) == ((1/12)*m*L^2)*theta_ddot;
eqn = transpose(eqn)

eqn = 

Non-Linear EOM

unknowns = solve(eqn,[y_ddot; theta_ddot]);
y_ddot = simplify(unknowns.y_ddot)

y_ddot = 

theta_ddot = simplify(unknowns.theta_ddot)

theta_ddot = 

Solving Nonlinear EOM

Make all state variables functions of time and declaring constants

L = 0.45;               % Length of the Pendulum
m = 1;                  % Mass of the Pendulum
g = 9.81;               % Gravitational Constant
y_0 = 0.005;            % Initial y
omega1 = 200*pi;        % 1st Omega in rad/s
omega2 = 140*pi;        % 2nd Omega in rad/s
theta_0 = 5*pi/180;     % Initial theta in rad
t_0 = 0;                % Initial time
t_f = 10;               % Final time
n = 300;                % Iterations

syms theta(t) y(t) theta_dot(t) y_dot(t) theta_ddot1(t) y_ddot1(t) theta_ddot2(t) y_ddot2(t)

Substitute these symbolic functions into the EOMs. We then have MATLAB produce the EOM function using the ode command.

theta_ddot1 = subs(theta_ddot, [str2sym('omega'), str2sym('y_0'), str2sym('m'), str2sym('L'), str2sym('g')], [omega1, y_0, m, L, g])

theta_ddot1 = 

theta_ddot2 = subs(theta_ddot, [str2sym('omega'), str2sym('y_0'), str2sym('m'), str2sym('L'), str2sym('g')], [omega2, y_0, m, L, g]);
y_ddot1 = subs(y_ddot, [str2sym('omega'), str2sym('y_0'), str2sym('m'), str2sym('g')], [omega2, y_0, m, g])

y_ddot1 = 

y_ddot2 = subs(y_ddot, [str2sym('omega'), str2sym('y_0'), str2sym('m'), str2sym('g')], [omega2, y_0, m, g]);

eom1 = odeFunction([theta_dot; y_dot; theta_ddot1; y_ddot1], [theta; y; theta_dot; y_dot], g, L, m, y_0); 

Error using mupadengine/feval2sym_NaNsingularity
Invalid parameter '981/100'.

Error in sym/odeFunction (line 118)
A = feval2sym_NaNsingularity(symengine, 'daetools::odeFunction', expr, vars, params{:});

%HERE is where I get the errors

eom2 = odeFunction([theta_dot; y_dot; theta_ddot2; y_ddot2], [theta; y; theta_dot; y_dot], g, L, m, y_0); 

[T1 S1] = ode45(@(t,s)eom1(t, s, L, m, y_0), linspace(t_0,t_f,n), [theta_0, y_0, 0, 0]);
[T2 S2] = ode45(@(t,s)eom2(t, s, L, m, y_0), linspace(t_0,t_f,n), [theta_0, y_0, 0, 0]);

Answer 1

[Steven Lord]

When you call odeFunction with g as part of the third input, it is not a symbolic variable. You originally defined it using syms but afterwards it was overwritten with the numeric value 9.81. Therefore it's not usable as a symbolic parameter of the system.