A self-consistent field (SCF) London-orbital computational scheme to perform gauge-origin independent nonperturbative calculations for molecules in strong magnetic fields is presented. The crucial difference in the proposed approach with respect to common-origin finite-field SCF implementations consists in the evaluation of molecular integrals over the field-dependent molecular basis functions, which is tantamount to computing molecular integrals in a hybrid Gaussian and plane-wave basis set. The implementation of a McMurchie-Davidson scheme for the calculation of the molecular integrals over London orbitals is discussed, and preliminary applications of the newly developed code to the calculation of fourth-rank hypermagnetizabilities for a set of small molecules, benzene, and cyclobutadiene are presented. The nonperturbative approach is particularly useful for studying the highly nonlinear response of paramagnetic closed-shell systems such as boron monohydride, or the pi-electron response of cyclobutadiene.