@abdgafartunde  Instead of column matrices, it is better to use vectors.

See the toy example:

N:= 5:
A:=<1,2; 3,4>;
b:=<1, 2>;  C:=<1, 1>;

for i from 0 to N do
      x[i+1] := A%.((i^2+1)*b) + C;
od; 

@abdgafartunde  I meant that in the same way, you can assign a name to each output, and then refer to these names when necessary.

@Ahmed111  See update to my answer above. The commands  diff  and  Physics:-diff  calculate any derivatives, both ordinary and partial. See help on these commands.

@Carl Love  See update to my answer. Compare my result with yours. 

@Scot Gould  Thank you for this!  I did not know that Maple solves systems written in vector form.

@amirhadiz  The  Interpolation  package appeared only in the latest versions of Maple. You have an older version, but you can use the code from my answer.

@acer  Thanks for this info. I did not know that such an evident command appeared recently only. Of course, for older versions, we can write
L1:=map(t->[t,ListTools:-Occurrences(t,L)], convert(L,set));
instead of  
L1:=ListTools:-Collect(L);
 but of course, Statistics:-Tally is more convenient. My point was simply to show that there are different ways to solve the same problem.

@Ali2020  It works in Maple 2019 only.

@minhthien2016 

Here is the equation with 6 another solutions:

abs(a*x+b)+abs(c*x+d)-3*x^2+n*x+p=0

@minhthien2016  
 

f:=x-> abs(a*x+b)+abs(c*x+d)-x^2+n*x+p;
solve([f(1) = 0, f(2) = 0, f(3) = 0, f(4) = 0, f(5) = 0, f(6) = 0], [a, b, c, d, n, p]);

     [[a = -2, b = 5, c = -2, d = 9, n = 7, p = -16], [a = -2, b = 5, c = 2, d = -9, n = 7, p = -16], [a = 2, b = -5, c = -2, d = 9, n = 7, p = -16], [a = 2, b = -5, c = 2, d = -9, n = 7, p = -16], [a = -2, b = 9, c = -2, d = 5, n = 7, p = -16], [a = 2, b = -9, c = -2, d = 5, n = 7, p = -16], [a = -2, b = 9, c = 2, d = -5, n = 7, p = -16], [a = 2, b = -9, c = 2, d = -5, n = 7, p = -16]]

@mehran rajabi  See help on  Student:-NumericalAnalysis:-Quadrature  command, spesifically the  method  option.

Present the complete code in an editable form (as text, not a picture) or the worksheet.

@weidade37211  I do not know how to explain this difference. I tried both methods and found that the second one is faster and does not need square brackets.

@acer  My comment above related to calculation of  cos(Pi/7)  and I used the de Moivre's formula and then  solve . The value of  sin(Pi/7)  probably cannot be simplified as well as  cos(Pi/7) .

@vv  I meant real radicals (without the imaginary unit  I ). Your expression can be slightly simplified, for example

cos(Pi/7) = (1/12)*(-28+(84*I)*sqrt(3))^(1/3)+7/(3*(-28+(84*I)*sqrt(3))^(1/3))+1/6

without  sqrt .