I have just found that (3) and (4) can be done by the following code
TRUE = 1;
FALSE = 0;
OK = FALSE;
while OK == FALSE
fprintf(1,'Input the number of equations\n');
M = input(' ');
if M <= 0 || M > 7
fprintf(1,'Number must be a positive integer < 8\n');
else
OK = TRUE;
end
end
ss = cell(M,1);
for I = 1:M
fprintf(1,'Input the function F_(%d) in terms of t and y1 ... y%d\n',...
I,M);
fprintf(1,'For example: ''y1-t^2+1'' \n');
kk = input(' ');
ss{I} = kk;
end
%(3)
V = 1:M+1; %replacing V with other input vectors here
c = num2cell(V);
Y = zeros(1,M);
if M==1
ysym = sym('y1');
else
ysym = sym('y%d',[1,M]);
end
syms t
for I = 1:M
z = evalin(symengine, ss{I});
f = symfun(z,[t,ysym]);
Y(I)=f(c{:});
end
%(4)
Y