The commercial macsyma has a function odematsys,
which I implemented, with the addition of accepting
initial conditions as a optional third argument.
I would appreciate comments and or tests of it.
Dan Stanger
load(desoln);
LOAD('LRATS);
odematsys(mateqn, vars, [intconds]):=block([indepvar,lvars,eqns,sol,j],
mateqn:lhs(mateqn)-rhs(mateqn),
indepvar:first(listofvars(vars)),
lvars:[],
for i in vars do lvars:append(lvars,i),
eqns:[],
for i in mateqn do eqns:append(eqns,i),
eqns:desolve(eqns,lvars),
sol:matrixmap(lambda([x],assoc(x,eqns)),vars),
j:0,
for intcond in intconds do block([lcond],
lcond:[],
for i in intcond do lcond:append(lcond,i),
sol:lratsubst(
maplist(lambda([x,y],at(diff(x,indepvar,j),indepvar=0)=y),lvars,lcond),
sol),
j:j+1
),
sol
);
m1:genmatrix(lambda([i,j],u[i,j](t)),2,2);
m2:genmatrix(lambda([i,j],'diff(u[i,j](t),t)),2,2);
m3:matrix([3,1],[1,3]);
eq:m2=m3.m1;
odematsys(eq,m1,matrix([1,0],[0,1]));