Recently we proposed a scheme that applies the variational principle to a coherent-pair condensate in the BCS case [1]. This work extends the scheme to the HFB case by allowing variation of the canonical single-particle basis. The result is equivalent to that of the so-called variation after particle-number projection in the HFB case, but now the particle number is always conserved and the time-consuming projection is avoided. Specifically, we derive the analytical expression for the gradient of the average energy with respect to changes of the canonical basis, which is then used in the family of gradient minimizers to iteratively minimize energy. In practice, we find the so-called ADAM minimizer (adaptive moments of gradient), borrowed from the machine-learning field, is very effective. The new algorithm is demonstrated in a semi-realistic example using the realistic V low-k interaction and large model spaces (up to 15 harmonic-oscillator major shells). It easily runs on a laptop, and practically the computer time cost (in the HFB case) is several times that of solving HF by the gradient minimizers. We hope the new algorithm could become a common practice because of its simplicity.