  Let $z=x+iy \in \mathbb{H}:=\{z= x+ i y\in\mathbb{C}: y>0\}$ and $ \theta
(\alpha;z)=\sum_{(m,n)\in\mathbb{Z}^2 } e^{-\alpha \frac{\pi }{y }|mz+n|^2}$ be
the theta function associated with the lattice $L ={\mathbb Z}\oplus z{\mathbb
Z}$. In this paper we consider the following minimization problem of difference
of two theta functions \begin{equation}\aligned\nonumber \min_{ \mathbb{H} }
\Big(\theta (\alpha; z)-\beta\theta (2\alpha; z)\Big) \endaligned\end{equation}
where $\alpha \geq 1$ and $ \beta \in (-\infty, +\infty)$. We prove that there
is a critical value $\beta_c=\sqrt2$ (independent of $\alpha$) such that if
$\beta\leq\beta_c$, the minimizer is $\frac{1}{2}+i\frac{\sqrt3}{2}$ (up to
translation and rotation) which corresponds to the hexagonal lattice, and if
$\beta>\beta_c$, the minimizer does not exist. Our result partially answers
some questions raised in \cite{Bet2016, Bet2018, Bet2020, Bet2019AMP} and gives
a new proof in the crystallization of hexagonal lattice under Yukawa potential.
