Needs["Graphics`PlotField`"]
<< graphics`ImplicitPlot`

Nonhomogeneous System

deq1  =  x1'[t]  == -2 x1[t]  +  x2[t]  +  
          Exp[-t]
deq2  =  x2'[t] == x1[t] - 2  x2[t] + 3  t

DSolve[{deq1, deq2}, {x1[t], x2[t]}, t]

Nonhomogeneous System with initial conditions

DSolve[{deq1, deq2, x1[0] == 0, x2[0] == 3}, {x1[t], 
        x2[t]}, t]

Equilibrium Solutions

Solve[{-2 x1[t]  +  x2[t] + Exp[-t] == 0, 
        x1[t] - 2  x2[t] + 3  t == 0}, {x1[t], x2[t]}]

steady   = 
      ParametricPlot[{(1/3) (E^(-t)) (4 + 
                3 (E^t) t), 
          (2/3) (E^(-t)) (1 + 
                3 (E^t) t)}, {t, 0, 2}, Frame -> True, 
        Axes -> False, 
        FrameLabel -> {"x1", "x2"}]

Nonlinear Homogeneous System

deq1  =  x1'[t]  == -5 x1[t]  + 2 x2[t]  +  
          x1[t] x2[t]
deq2  =  x2'[t] == 
        x1[t] + x2[t] - 4  x2[t] x1[t]

DSolve[{deq1, deq2}, {x1[t], x2[t]}, t]

Mathematica was unable to solve this analytically. But we can still learn something about the solution, and could always get a numerical solution if we wished.

Equilibrium Solutions

Clear[f1, f2]

f1[x1_, x2_]  =  -5 x1  + 2 x2 +  
        x1 x2
f2[x1_, x2_]  =  x1 +  x2  - 4  x2 x1

Solve[{f1[x1, x2] == 0 , 
        f2[x1, x2] == 0}, {x1, x2}]

steady  =  
      ImplicitPlot[{f1[x1, x2] == 0 , 
          f2[x1, x2] == 0}, {x1, -10, 10}, {x2, -2, 10}, 
        Frame -> True, FrameLabel -> {"x1", "x2"}, 
        PlotStyle -> { Dashing[{1.0, 0.0}], 
            Dashing[{0.01, 0.01}]}]

field = 
      PlotVectorField[{f1[x1, x2], f2[x1, x2]}, {x1, -1, 1}, {x2, -1, 
          1}, Frame -> True, 
        FrameLabel -> {"x1", "x2"}]
Show[field, Graphics[Disk[{0, 0}, 0.05]], 
      Graphics[Disk[{7/19, 7/9}, 0.05]]]

field = 
      PlotVectorField[{f1[x1, x2], f2[x1, x2]}, {x1, -15, 
          15}, {x2, -10, 10}, Frame -> True, 
        FrameLabel -> {"x1", "x2"}, 
        PlotPoints -> 6]
steady  =  
      ImplicitPlot[{f1[x1, x2] == 0 , 
          f2[x1, x2] == 0}, {x1, -15, 15}, {x2, -10, 10}, 
        Frame -> True, FrameLabel -> {"x1", "x2"}, 
        PlotStyle -> { Dashing[{1.0, 0.0}], 
            Dashing[{0.01, 0.01}]}]
Show[field, Graphics[Disk[{0, 0}, 0.05]], 
      Graphics[Disk[{7/19, 7/9}, 0.05]], steady]

field = 
      PlotVectorField[{f1[x1, x2], f2[x1, x2]}, {x1, -1/4, 
          1/4}, {x2, -1/4, 1/4}, Frame -> True, 
        FrameLabel -> {"x1", "x2"}]
Show[field, Graphics[Disk[{0, 0}, 0.02]]]
field = 
      PlotVectorField[{f1[x1, x2], f2[x1, x2]}, {x1, 0, 1/2}, {x2, 1/2, 1}, 
        Frame -> True, 
        FrameLabel -> {"x1", "x2"}]
Show[field, Graphics[Disk[{7/19, 7/9}, 0.02]]]

Sensitivity Analysis (Ignored Terms)

f1[x1_, x2_]  =  -5 x1  + 2 x2 +  x1 x2  -  
        10^(-1) x2^2
f2[x1_, x2_]  =  x1 +  x2  - 4  x2 x1
Solve[{f1[x1, x2] == 0 , f2[x1, x2] == 0}, {x1, 
        x2}]
N[%]

field = 
      PlotVectorField[{f1[x1, x2], f2[x1, x2]}, {x1, 0, 0.75}, {x2, 21.5, 
          22.25}, Frame -> True, 
        FrameLabel -> {"x1", "x2"}, 
        PlotPoints -> 10]
steady  =  
      ImplicitPlot[{f1[x1, x2] == 0 , f2[x1, x2] == 0}, {x1, 0, 
          1}, {x2, 21, 23}, Frame -> True, 
        FrameLabel -> {"x1", "x2"}, 
        PlotStyle -> { Dashing[{1.0, 0.0}], 
            Dashing[{0.01, 0.01}]}]
Show[field, Graphics[Disk[{0.25, 21.95}, 0.02]]]