LCAO is essentially a variational method, where all we are doing is guessing the answer and calculating the energies that result. These energies may or may not be close to the actual experimental answers, which we should always refer to.
Let me briefly derive the answer for $\text{H}2^+$ using the full method... As you say, we start by saying we would like a solution of the form $\psi=c_1\psi\text{1s}^a+c_2\psi_\text{1s}^b$. 976
In the variational method, we say that the quantity given by $$\widetilde{E}=\frac{\int\psi^\hat{H}\psi\;\text{d}\tau}{\int\psi^\psi\;\text{d}\tau}$$ is the energy of our wavefunction. $\hat{H}$ is the Hamiltonian for the system, and $\text{d}\tau$ indicates we should integrate over all space. Sticking in $\psi$ and after a bit of rearranging, we get: $$\widetilde{E}=\frac{c_1^2\int\psi_\text{1s}^a\hat{H}\psi_\text{1s}^a\;\text{d}\tau+c_2^2\int\psi_\text{1s}^b\hat{H}\psi_\text{1s}^b\;\text{d}\tau+2c_1c_2\int\psi_\text{1s}^a\hat{H}\psi_\text{1s}^b\;\text{d}\tau}{c_1^2\int\psi_\text{1s}^a\psi_\text{1s}^a\;\text{d}\tau+c_2^2\int\psi_\text{1s}^b\psi_\text{1s}^b\;\text{d}\tau+2c_1c_2\int\psi_\text{1s}^a\psi_\text{1s}^b\;\text{d}\tau}$$ Note I've assumed that everything is real, and that for example $\int\psi_\text{1s}^b\hat{H}\psi_\text{1s}^a\;\text{d}\tau=\int\psi_\text{1s}^a\hat{H}\psi_\text{1s}^b\;\text{d}\tau$. Let's define some constants to make the integrals go away: \begin{align} \alpha&=\int\psi_\text{1s}^a\hat{H}\psi_\text{1s}^a\;\text{d}\tau=\int\psi_\text{1s}^b\hat{H}\psi_\text{1s}^b\;\text{d}\tau\ \beta&=\int\psi_\text{1s}^b\hat{H}\psi_\text{1s}^a\;\text{d}\tau=\int\psi_\text{1s}^a\hat{H}\psi_\text{1s}^b\;\text{d}\tau\ S&=\int\psi_\text{1s}^a\psi_\text{1s}^b\;\text{d}\tau \end{align} These can be interpreted as the self-energy ($\alpha$), the interaction energy ($\beta$) and the orbital overlap ($S$). 526