Weird Matlab error with piecewise element in matrix, changes matrix dimension

I'm trying to solve a differential equation that has the form Y'(t)=A(t)*Y(t), where Y and Y' is a column with 4 elements, A is a 4x4 matrix. In the A(t) matrix I'm using the functions m(t) and f(t) to compute each value of the matrix. However, when I put a simple value to m(t) and f(t) (ex: a constant value or a t^2...) the function solvedPmWithTime returns a 4x1 matrix, which is what I expect, but when I use the m(t) and /or function f(t) that I want, which is computed by the compute_m_and_f function, it returns a 1x1 matrix. I can't find my error. If someone could help it would be very much appreciated :). Here are the codes of the functions and the code of the live script that I use to compute the differential equation:

b normally is 4x4 but as you see it sais 1x1

function [m,f] = compute_m_and_f()
syms m(t) f(t)
Bm =600000;
sigmam =2.5;
tm =6;
tf =6;
sigmaf = 2;
Bf = 600000;

eqnM = diff(m,t) == Bm / sqrt(2*pi*sigmam^2)*exp(-((t-tm)^2)/(2*sigmam^2));
condM = m(0)==0;

eqnF = diff(f,t) == Bf / sqrt(2*pi*sigmaf^2)*exp(-((t-tf)^2)/(2*sigmaf^2));
condF = f(0)==0;

M(t) = dsolve(eqnM, condM);
F(t) = dsolve(eqnF, condF);

taum = 4;
tauf = 5;

m(t) = piecewise(t>taum, M(t)-M(t-taum),M(t));
f(t) = piecewise(t>tauf, F(t)-F(t-taum),F(t));
end

function solvedPmWithTime= solvedPmWithTime(~, y)

[m, f] = compute_m_and_f;

alpha = 480;
beta = 0.0011;

sigma = 0.0634;
epsilon = 0.0725;

a = (alpha * (f + m)) / (1 + alpha * beta * (f + m));

b11 =- (1 - (m / (f + m)) ^ sigma);
b12 = (1 - (f / (f + m)) ^ epsilon);
b13 = (1 - (f / (f + m)) ^ epsilon);
b14 = (1 - (f / (f + m)) ^ epsilon);

b21 = (1 - (m / (f + m)) ^ sigma);
b22 = -1;
b23 = 0;
b24 = 0;

b31 = 0;
b32 = (f / (f + m)) ^ epsilon;
b33 = -1;
b34 = 0;

b41 = 0;
b42 = 0;
b43 = (f / (f + m)) ^ epsilon;
b44 =- (1 - (f / (f + m)) ^ epsilon);

b = [b11, b12, b13, b14; b21, b22, b23, b24; b31, b32, b33, b34; b41, b42, b43, b44];
A = a * b;
solvedPmWithTime = A * y;
end

clear all
rho = 6;
Bf = 600000;

Pm1_0 = 0;
Pm2_0 = 0;
Pm3_0 = 0;
Pf_0 = 3600;
[t, Pm] = ode45(@solvedPmWithTime, [0, 100], [Pm1_0; Pm2_0; Pm3_0; Pf_0]); %this returns an error because the size of the matrix is not 4x1

Pm1 = Pm(:, 1);
Pm2 = Pm(:, 2);
Pm3 = Pm(:, 3);
Pf = Pm(:, 4);
hold off
plot(t, Pm1, t, Pm2, t, Pm3, t, Pf)
legend('Pm1', 'Pm2', 'Pm3', 'Pf')



n = solvedPmWithTime([0, 100], [Pm1_0; Pm2_0; Pm3_0; Pf_0]);
size(n)

Answer 1

[Rasmita]

Hi,

It is my understanding that, you are trying to solve a differential equation involving matrices. The function ‘solvedPmWithTime’ should ideally return a 4x1 matrix, but it’s returning a 1x1 matrix.

The reason for the error is that the 'compute_m_and_f' function returns symbolic functions ‘m(t)’ and ‘f(t)’ instead of their numerical values. Therefore, when you use these functions in the ‘solvedPmWithTime’ function, the resulting matrix ‘A’ contains symbolic values, which causes the ‘ode45’ solver to fail.

To fix this, you need to evaluate the symbolic functions ‘m(t)’ and ‘f(t)’ at the corresponding time points in the ‘solvedPmWithTime’ function. You can do this using the ‘subs’ function, which substitutes numerical values for symbolic variables. So, replace ‘m’ and ‘f’ in ‘solvedPmWithTime’ function ‘a’ and ‘b’ matrix elements calculation with ‘subs(m)’ and ‘subs(f)’.

Also use ‘double’ function to convert the resulting symbolic matrix ‘A’ to a numerical matrix in function ‘solvedPmWithTime’, as follows:

A = double(a * b);

For more information, please refer below documentation links:

https://www.mathworks.com/help/symbolic/subs.html

https://www.mathworks.com/help/matlab/ref/double.html

Hope this helps!

Regards,

Rasmita