We present and analyze a uniquely solvable and unconditionally energy stable numerical scheme for the Functionalized Cahn–Hilliard equation, including an analysis of convergence. One key difficulty associated with the energy stability is based on the fact that one nonlinear energy functional term in the expansion is neither convex nor concave. To overcome this subtle difficulty, we add two auxiliary terms to make the combined term convex, which in turns yields a convex–concave decomposition of the physical energy. As a result, both the unconditional unique solvability and the unconditional energy stability of the proposed numerical scheme are assured. In addition, a global in time $$H_{\mathrm{per}}^2$$ Hper2 stability of the numerical scheme is established at a theoretical level, which in turn ensures the full order convergence analysis of the scheme, which is the first such result in this field. To deal with an implicit 4-Laplacian term at each time step, we apply an efficient preconditioned steepest descent algorithm to solve the corresponding nonlinear systems in the finite difference set-up. A few numerical results are presented, which confirm the stability and accuracy of the proposed numerical scheme.