Many important problems in Operations Research and Statistics require the computation of nondominated (or Pareto or efficient) sets. This task may be currently undertaken efficiently for discrete sets of alternatives or for continuous sets under special and fairly tight structural conditions. Under more general continuous settings, parametric characterisations of the nondominated set, for example through convex combinations of the objective functions or -constrained problems, or discretizations-based approaches, pose several problems. In this paper, the lack of a general approach to approximate the nondominated set in continuous multiobjective problems is addressed. Our simulation-based procedure only requires to sample from the set of alternatives and check whether an alternative dominates another. Stopping rules, efficient sampling schemes, and procedures to check for dominance are proposed. A continuous approximation to the nondominated set is obtained by fitting a surface through the points of a discrete approximation, using a local (robust) regression method. Other actions like clustering and projecting points onto the frontier are required in nonconvex feasible regions and nonconnected Pareto sets. In a sense, our method may be seen as an evolutionary algorithm with a variable population size.