solving system of ODE and exporting results into matrix form

Hi, all.

I am new with Maxima and I'd like to use it for automation of the
following job:

I have a system of ODE describing the simplest differential-difference
heat-transfer model in sensor:
eq1:q1(tau)=c*ro*S*delta/2*diff(t1(tau),tau)+lambda*S*(t1(tau)-t2(tau))/delta;
eq2:lambda*S*(t1(tau)-t2(tau))/delta=c*ro*S*delta*diff(t2(tau),tau)+lambda*S*(t2(tau)-t3(tau))/delta;
eq3:lambda*S*(t2(tau)-t3(tau))/delta=c*ro*S*delta/2*diff(t3(tau),tau)+q2(tau);

I can solve equations relative to derivative:
s1:solve(eq1, diff(t1(tau),tau));
s2:solve(eq2, diff(t2(tau),tau));
s3:solve(eq3, diff(t3(tau),tau));

But I can't do following:

1) format results in form "diff(t1(tau),tau)
=A*t1(tau)+B*t2(tau)+C*q1(tau)"

2) substitute variables in results "lambda/(c*ro*delta^2)" with "b" and
"1/(c*ro*delta)" with "d"

3) represent system of ODE in "dt_dtau=F*t+G*u" matrix form, where:
dt_dtau:transpose(matrix([diff(t1(tau),tau),diff(t2(tau),tau),diff(t3(tau),tau)]));
t:transpose(matrix([t1(tau),t2(tau),t3(tau)]));
u:transpose(matrix([q1(tau),q2(tau)]));

Thanks in advance.

-- 
Kirill Kirillov
<><