@Markiyan Hirnyk 

Thanks for the comment. Indication of editing done.

@oldstudent 

with(Student[Calculus1]):

A := diff(x^2*sin(x), x);

ShowSolution(A);

Another way to solve the problem, first create a matrix of the correct size, and then make the necessary substitutions or conversions.

x := [i $ i=1..25];

A:=Matrix(5, x);

eval(A, [seq(x[i]=convert(x[i], string), i=1..nops(x))]);

@bluehotel  Just copy the bold text below into Mathematica and press the key Enter:

Reduce [{-x3 * x5 /m - x3 * x5 - 2 * x1 * x5 + 2 * x1 == 0, 
x3 * x5 /m - 2 * x2 * x5 + 2 * x2 == 
0,-х3 + 2 * x5 - x3 * x5 + 1 - 3 * x1 - 3 * x2 + x4 == 
0, -2 * x4 - 2 * x4 * x5 == 0,-x3 * x5 ^ 2 / (x1 * m) - 4 * х5 ^ 2 + 4 * x5 == 
0}, {x1, x2, x3, x4, x5}] // TraditionalForm

@liumai   rhs((P(k))(t)[2])  is the right side of the second operand of the list  (P(k))(t)  for specific  k  and  t .

For 2 functions in the system replace  (k, t)->rhs((P(k))(t)[2])  by the list  [(k, t)->rhs((P(k))(t)[2]),(k, t)->rhs((P(k))(t)[4])]

@Carl Love  Your code in M 12 also works. I think  evala  command is crucial. I used  expand  command instead  and got the message "lost connection with the kernel"

@Alejandro Jakubi

1)  Thanks for the useful information regarding  parametric. 

2) You replaced the equation  sqrt(x-a) = x  with  the equation  x-a = x^2.  But  in the real domain, these equations are not equivalent.

@Carl Love  In my code the specific function  f  was used rather than arbitrary function as in your code. I think the questioner himself according this example will write the function that interests him.

@mahmood180  Your example for  tj=1:

HybrFunc(4, 3, 1):

for j from 0 to 2 do

seq('b'[i,j]=b[i,j](t), i=1..4);

od;

Download the text (not an image) of the code, so we can check its accuracy.

@Stephan   

1) See the help on  piecewise

2) Add the line

applyop(x -> lhs(x)/L < rhs(x)/L, {1,3}, Mf(xi));

Animation quality can be slightly improved:

PDE := diff(u(x, t), t$2) = diff(u(x, t), x$2):

IBCs := u(x, 0) = sin((1/2)*Pi*x)*exp(-x), (D[2](u))(x, 0) = 0, u(0, t) = 0, u(4, t) = 0:

Sol := pdsolve({PDE}, {IBCs}, numeric):

Sol:-animate(t = 0 .. 20, frames = 300, thickness = 2, numpoints = 3000, labels=[x, u(x, t)]);

@Markiyan Hirnyk  There are infinitely many matrices  A  and  B  satisfying the condition  A.B = C . This is seen from the direct solution:

restart;

C:=Matrix([[8, 2, -2], [2, 5, 4], [-2, 4, 5]]):

A:=Matrix(3,2, symbol=a):

B:=Matrix(2,3, symbol=b):

solve(Equate(A.B,C)):

assign(%):

'A'=A; 'B'=B;

simplify(B.A);

@Joe Riel  of siderals is  Equate(A.B, C)

@Markiyan Hirnyk   I think that this is a natural way to solve such problems.