Why do we chose exponential function as a trial solution for second order linear differential equation with constant coefficient?

Consider a homogeneous linear second order ODE with constant coefficients, say $$ay''+by'+cy=0.$$ Denote by $E$ the vector space of real-valued functions of class $C^\infty$ on $\mathbb{R}$ and define the endomorphism $\partial$ of $E$ as $$\forall f\in E,\ \partial f=f'$$ (yes, just the derivative).

Our equation can be written as $$\bigl(a\partial^2+b\partial+c\text{id}\bigr)y=0,$$ (where $\text{id}$ denotes the identity endomorphism of $E$) so the solutions can be seen as the kernel of the endomorphism $\bigl(a\partial^2+b\partial+c\bigr)$.

Now factor (in $\mathbb{C}$ if needed) the polynomial $aX^2+bX+c$, say $$aX^2+bX+c=a(X-\alpha_1)(X-\alpha_2).$$ It is straightforward to check that $$a\partial^2+b\partial+c\text{id}=a(\partial-\alpha_1\text{id})\circ(\partial-\alpha_2\text{id}).$$ If you hit non-real numbers in the previous factorization, you might want to extend your field of scalars to $\mathbb{C}$ and consider $E\otimes_{\mathbb{R}}\mathbb{C}$ instead of $E$.

The eigenvectors of $\partial$ being the exponentials, it's quite natural to use them.

In fact, this factorization of the pseudo-differential operator in the spirit of Weyl algebras leads to a nice proof of the form of the solutions of the ODE.

With this point of view, it seems natural that there are no sensible trial functions other than the elements in the characteristic spaces of $\partial$, namely the functions of the form $$x\mapsto P(x)\mathrm{e}^{\lambda x}$$ where $P$ is a polynomial function and $\lambda\in\mathbb{R}$ or $\lambda\in\mathbb{C}$. In the case $\alpha_1=\alpha_2$ we actually see such solutions.

All this can be generalized to ODEs of any order, with constant coefficients; for non-constant coefficients you can, depending on how the coefficients are, do a similar study. With non-constant coefficients, the problem will be in factorizing the pseudo-differential operator.

The exponential function $e^{\lambda x}$ is special because its derivative is a multiple of itself. Thus, substituting it into the ODE gives a quadratic equation in $\lambda$, which is easy to solve. In general, this method gives us ALL the possible solutions to the ODE, and so using other functions as guess solutions won't work.

If $y'=ay$ consider $z=e^{-ax}y$. We have $z'=-ae^{-ax}y+e^{-ax}y'=e^{-ax}(y'-ay)=0$

It is elementary then that $z=A$ where $A$ is a constant and therefore $y=Ae^{ax}$

Now suppose $$y''-(a+b)y'+aby=0$$ and take $p=y'-ay; q=y'-by$ so that we find $$p'-bp=0; q'-aq=0$$

whence (from the first result) $$y'-ay=p=Ae^{bx}; y'-by=q=Be^{ax}$$

Now considering $p-q$ we eliminate $y'$ to give $$(b-a)y=Be^{ax}-Ae^{bx}$$Which is a solution of the form you wanted - and the logic (which depends on the mean value theorem) shows that this is the inevitable form of the solution to the equation, which is why we choose it as a test solution and determine the constants from the boundary conditions if necessary.

Obviously this goes wrong when $a=b$, but I will leave you to explore that case yourself.