Gridless methods show great superiority in line spectral estimation. These methods need to solve an atomic $l_0$ norm (i.e., the continuous analog of $l_0$ norm) minimization problem to estimate frequencies and model order. Since this problem is NP-hard to compute, relaxations of atomic $l_0$ norm, such as nuclear norm and reweighted atomic norm, have been employed for promoting sparsity. However, the relaxations give rise to a resolution limit, subsequently leading to biased model order and convergence error. To overcome the above shortcomings of relaxation, we propose a novel idea of simultaneously estimating the frequencies and model order by means of the atomic $l_0$ norm. To accomplish this idea, we build a multiobjective optimization model. The measurment error and the atomic $l_0$ norm are taken as the two optimization objectives. The proposed model directly exploits the model order via the atomic $l_0$ norm, thus breaking the resolution limit. We further design a variable-length evolutionary algorithm to solve the proposed model, which includes two innovations. One is a variable-length coding and search strategy. It flexibly codes and interactively searches diverse solutions with different model orders. These solutions act as steppingstones that help fully exploring the variable and open-ended frequency search space and provide extensive potentials towards the optima. Another innovation is a model order pruning mechanism, which heuristically prunes less contributive frequencies within the solutions, thus significantly enhancing convergence and diversity. Simulation results confirm the superiority of our approach in both frequency estimation and model order selection.