reset();

solucoes_fundamentais(a,b,c,t):=block
(
  [delta,y1,y1,alpha,beta],
  delta:b^2-4*a*c,
  if delta>0 then
    (
      y1:exp((-b+sqrt(delta))*t/(2*a)),
      y2:exp((-b-sqrt(delta))*t/(2*a))
    )
  elseif delta=0 then
    (
      y1:exp(-b*t/(2*a)),
      y2:t*exp(-b*t/(2*a))
    )
  else
    (
      alpha:-b/(2*a),
      beta:sqrt(-delta)/(2*a),
      y1:cos(beta*t)*exp(alpha*t),
      y2:sin(beta*t)*exp(alpha*t)
    ),
  [y1,y2]
);

y:solucoes_fundamentais(1,5,6,x)$display(y[1],y[2]);

y:solucoes_fundamentais(1,2,1,x)$display(y[1],y[2]);

y:solucoes_fundamentais(1,2,5,x)$display(y[1],y[2]);

wrosk(y1,y2,t):=block
(
  [M],
  M:matrix([y1,y2],[diff(y1,t),diff(y2,t)]),
  determinant (M)
);

W:wrosk(y[1],y[2],x);

trigsimp(W);

yp(y,f,t):=block
(
  [W,int1,int2],
  W:wrosk(y[1],y[2],t),
  int1:integrate(y[2]*f/W,t),
  int2:integrate(y[1]*f/W,t),
  -y[1]*int1+y[2]*int2
);

y:solucoes_fundamentais(1,5,6,x)$display(y[1],y[2]);

sol_par:yp(y,x*exp(-5*x),x);

ratsimp(sol_par);

geral_2_ordem(a,b,c,f,t,C1,C2):=block
(
  [y,YP],
  y:solucoes_fundamentais(a,b,c,t),
  YP:yp(y,f,t),
  C1*y[1]+C2*y[2]+YP
);

geral_2_ordem(1,5,6,x*exp(-5*x),x,C1,C2);

Created with wxMaxima.