restart:
with(linalg):
A:=matrix(2,2,[12,-2,-7,-1]);
eigenvalues(A);
eigenvectors(A);
Q:=matrix(2,2,[-2,1,1,7]);
Diag:=matrix(2,2,[13,0,0,-2]);
evalm(A&*Q);
evalm(Q&*Diag);
Qinv:=inverse(Q);
evalm(Q&*Diag&*Qinv);
Q2:=matrix(2,2,[1,2,-1/2,14]);
Q2inv:=inverse(Q2);
evalm(Q2&*Diag&*Q2inv);
Diag2:=matrix(2,2,[-2,0,0,13]);
QSwap:=matrix(2,2,[1,-2,7,1]);
QSwapInv:=inverse(QSwap);
evalm(QSwap&*Diag2&*QSwapInv);
eDt:=matrix(2,2,[exp(13*t),0,0,exp(-2*t)]);
Soln1:=collect(simplify(evalm(Q&*eDt&*Qinv&*matrix(2,1,[x0,y0]))),{exp(13*t),exp(-2*t)});
Xsoln:=Soln1[1,1];
Ysoln:=Soln1[2,1];
diff(Xsoln,t)-(12*Xsoln-2*Ysoln);
simplify(%);
diff(Ysoln,t)-(-7*Xsoln-1*Ysoln);
simplify(%);
F:=matrix(2,1,[sin(3*t),t*exp(-7*t)]);
emDt:=matrix(2,2,[exp(-13*t),0,0,exp(2*t)]);
G1:=evalm(emDt&*Qinv&*F)[1,1];
G2:=evalm(emDt&*Qinv&*F)[2,1];
newG1:=combine(expand(G1));
newG2:=combine(expand(G2));
int(newG1,t);
int(newG2,t);
intG:=matrix(2,1,[int(newG1,t),int(newG2,t)]);
C:=matrix(2,1,[Xc,Yc]);
Soln:=evalm(Q&*eDt&*intG+Q&*eDt&*Qinv&*C);
X:=simplify(Soln[1,1]);
Y:=simplify(Soln[2,1]);
simplify(diff(X,t)-(12*X-2*Y+sin(3*t)));
simplify(diff(Y,t)-(-7*X-1*Y+t*exp(-7*t)));
simplify(subs({t=0},[X,Y]));