  Within the framework of non-relativistic QCD (NRQCD) factorization, we
compute the matching coefficients between full Quantum Chromodynamics (QCD) and
NRQCD for the heavy flavor-changing vector, axial-vector, scalar and
pseudo-scalar currents up to next-to-next-to-next-to-leading order (N$^3$LO).
We accomplish the analytical expressions for the three-loop renormalization
constants and the corresponding anomalous dimensions for all of the four NRQCD
currents with two different heavy flavors. The three-loop QCD corrections to
the matching coefficients turn out to be significantly larger than lower order
corrections. By employing the scale relation, we obtain the N$^3$LO corrections
to the wave functions at the origin for the vector $B_c^*$ meson and the
pseudo-scalar $B_c$ meson from the known result for the equal-mass heavy
quarkonium in potential NRQCD (pNRQCD). We find large cancellations at the
third order between the matching coefficients and the wave functions at the
origin, and obtain the convergent decay constants of $B^*_{c}$ and $B_{c}$ up
to N$^3$LO. We present the complete perturbative QCD predictions for the decay
constants, leptonic decay widths, and branching ratios of the beauty-charmed
mesons.
