LCAO (Linear Combination of Atomic Orbitals) and Phases
The concept of linear combination of atomic orbitals (LCAO) to form molecular orbitals (MO) is probably best understood, while digging a little deeper into quantum chemistry. The method is an approximation that was introduced for ab initio methods like Hartree Fock. I don not want to go into too much detail, but there a some points that need to be considered before understanding what LCAO actually does.

The time independent Schrödinger equation can was postulated as $$\mathbf{H}\Psi=E\Psi,$$
 with the Hamilton operator $\mathbf{H}$, the wave function $\Psi$ and the corresponding energy eigenvalue(s) $E$. I will just note that we are working in the framework of the Born-Oppenheimer approximation and refer to many textbooks for more details.

There is a set of rules, the wave function has to obey. 


It is a scalar, that can be real or complex, but the product of itself with its complex conjugated version is always positive and real. $$0\leq \Psi^*\Psi=|\Psi|^2$$
The probability of finding all $N$ electrons in all space $\mathbb{V}$ is one, hence to function is normed. $$N=\int_\mathbb{V} |\Psi(\mathbf{x}_1,\mathbf{x}_2,\dots,\mathbf{x}_N)|^2 \mathrm{d}(\mathbf{x}_1,\mathbf{x}_2,\dots,\mathbf{x}_N)$$ 
The value of the wave function has to vanish at infinity. $$0 = \lim_{\mathbf{x}\to\infty}|\Psi(\mathbf{x})|$$
The wave function has to be continuous and continuously differentiable, due to the second order differential operator for the kinetic energy $\mathbf{T}_c$ included in $\mathbf{H}$.
The Pauli Principle has to be obeyed. $$\Psi(\mathbf{x}_1,\mathbf{x}_2,\dots,\mathbf{x}_N) = -\Psi(\mathbf{x}_2,\mathbf{x}_1,\dots,\mathbf{x}_N)$$


Also the variational principle should hold for this approximation, stating that the expectation value for the energy of any trial wave function is larger that the energy eigenvalue of the true ground state. One of the most basic methods to approximately solve this problem is Hartree Fock. In it the trial wave function $\Phi$ is set up as a Slater determinant.
$$
\Phi(\mathbf{x}_1, \mathbf{x}_2, \ldots, \mathbf{x}_N) =
\frac{1}{\sqrt{N!}}
\left|
   \begin{matrix} 
      \phi_1(\mathbf{x}_1) &amp; \phi_2(\mathbf{x}_1) &amp; \cdots &amp; \phi_N(\mathbf{x}_1) \\
      \phi_1(\mathbf{x}_2) &amp; \phi_2(\mathbf{x}_2) &amp; \cdots &amp; \phi_N(\mathbf{x}_2) \\
      \vdots               &amp; \vdots               &amp; \ddots &amp; \vdots \\
      \phi_1(\mathbf{x}_N) &amp; \phi_2(\mathbf{x}_N) &amp; \cdots &amp; \phi_N(\mathbf{x}_N)
   \end{matrix} \right| 
$$

The expectation value of the energy in the Hartree Fock formalism is set up, using bra-ket notation as $$E =\langle\Phi|\mathbf{H}|\Phi\rangle.$$
Skipping through some major parts of the deviation of HF, well end up at an expression for the energy.
$$E=\sum_i^N \langle\phi_i|\mathbf{H}^c|\phi_i\rangle +\frac12 \sum_i^N\sum_j^N \langle\phi_i|\mathbf{J}_j-\mathbf{K}_j|\phi_i\rangle$$

To find the best one-electron wave functions $\phi_i$ we introduce Lagrange multiplicators $\lambda$ minimising the energy with respect to our chosen conditions. These conditions include that the molecular orbitals are ortho normal. $$\langle\phi_i|\phi_j\rangle =\delta_{ij} =\left\{
\begin{matrix}  0 &amp; , \text{for}~i \neq j\\ 1 &amp; , \text{for}~i = j\\\end{matrix}\right.$$
I will again skip through most of it and just show you the end expression.
$$\sum_j\lambda_{ij}|\phi_j\rangle = \mathbf{F}_i|\phi_i\rangle$$
with the Fock operator set up as
$$\mathbf{F}_i = \mathbf{H}^c +\sum_j (\mathbf{J}_j − \mathbf{K}_j)$$
with $i\in1\cdots{}N$, the total number of electrons.
We can transform these trial wavefunctions $\phi_i$ to canonical orbitals $\phi_i'$ (molecular orbitals) and obtain the pseudo eigenwertproblem $$\varepsilon_i\phi_i'=\mathbf{F}_i\phi_i'.\tag{Fock}$$ 
This equation is actually only well defined for occupied orbitals and these are the orbitals that give the lowest energy. In practice this formalism can be extended to include virtual (unoccupied) molecular orbitals as well.

Until now we did not use any atomic orbitals at all. This is the next step to find an approximation to actually solve these still pretty complicated systems.
LCAO a superposition method. In this approach we map a finite set of $k$ atomic (spin) orbitals $\chi_a$ onto another finite set of $l$ molecular (spin) orbitals $\phi_i'$. They are related towards each other via the expression 
\begin{align}
\phi_i'(\mathbf{x}) &amp;= c_{1,i}\chi_1(\mathbf{x}) + c_{2,i}\chi_2(\mathbf{x}) + \cdots + c_{k,i}\chi_k(\mathbf{x})\\
\phi_i'(\mathbf{x}) &amp;= \sum_a^k c_{a,i}\chi_a(\mathbf{x})\\\tag{LCAO}
\end{align}

From $\text{(Fock)}$ you can see, that there will be $l$ possible equations depending on the chosen set of orbitals, in the form of
$$\sum_a^k\sum_b^k c_{ia}^* c_{ib} \langle\chi_a|\mathbf{F}_i|\chi_b\rangle = \varepsilon_i \sum_a^k\sum_b^k c_{ia}^* c_{ib} \langle\chi_a|\chi_b\rangle$$
or for short there are $l$ equations $$C_{ab}F_{ab}=\varepsilon_i C_{ab}S_{ab}.$$ Since $a,b \in [1,k]$ you can gather $C_{ab}$ in a $k\times k$ matrix $\mathbb{C}$. Since $i\in[1,l]$, there will only be $l$ $\varepsilon_i$ and therefore the matrix of the Fock elements $F_{ab}$ has to be a $l\times l$ matrix $\mathbb{F}$. The whole problem reduces a matrix equation $$\mathbb{F}\mathbb{C}=\mathbb{S}\mathbb{C}\epsilon\!\!\varepsilon,\tag{work}$$
from which it is obvious, that the dimension of the involved matrices have to be the same. Hence $k=l$.

Too long, didn't read

The total number of elements in the finite set of atomic orbitals is equal to the total number of elements in the finite set of molecular orbitals. Any linear combination is possible, but only the orbitals that minimise the energy will be occupied.

The coefficients $c$ from $\text{(LCAO)}$ will be chosen to minimise the energy of the wave function. This will always be constructive interference. This is also independent of the "original" phase of the orbital that is combined with another orbital. In other words, from $\text{(LCAO)}$
$$\phi_i'(\mathbf{x}) = c_{1,i}\chi_1(\mathbf{x}) + c_{2,i}\chi_2(\mathbf{x}) \equiv c_{1,i}\chi_1(\mathbf{x}) - c_{2,i}[-\chi_2(\mathbf{x})].$$