Context.High-resolution spectroscopy has revealed large concentrations of CNO and sometimes other intermediate-mass elements (e.g., Ne, Na, Mg, or Al, for ONe novae) in the shells ejected during nova outbursts, suggesting that the solar composition material transferred from the secondary mixes with the outermost layers of the underlying white dwarf during thermonuclear runaway.Aims.Multidimensional simulations have shown that Kelvin-Helmholtz instabilities provide self-enrichment of the accreted envelope with material from the outermost layers of the white dwarf, at levels that agree with observations. However, the Eulerian and time-explicit nature of most multidimensional codes used to date and the overwhelming computational load have limited their applicability, and no multidimensional simulation has been conducted for a full nova cycle.Methods.This paper explores a new methodology that combines 1D and 3D simulations. The early stages of the explosion (i.e., mass-accretion and initiation of the runaway) were computed with the 1D hydrodynamic codeSHIVA. When convection extended throughout the entire envelope, the structures for each model were mapped into 3D Cartesian grids and were subsequently followed with the multidimensional codeFLASH. Two key physical quantities were extracted from the 3D simulations and were subsequently implemented intoSHIVA, which was used to complete the simulation through the late expansion and ejection stages: the time-dependent amount of mass dredged-up from the outer white dwarf layers, and the time-dependent convective velocity profile throughout the envelope.Results.This work explores for the first time the effect of the inverse energy cascade that characterizes turbulent convection in nova outbursts. More massive envelopes have been found that are those reported from previous models with pre-enrichment. These result in more violent outbursts, characterized by higher peak temperatures and greater ejected masses, with metallicity enhancements in agreement with observations.