We present an algorithm for constructing efficient surrogate frequency-domain models of (nonlinear) parametric dynamical systems in a non-intrusive way. To capture the dependence of the underlying system on frequency and parameters, our proposed approach combines rational approximation and smooth interpolation. In the approximation effort, locally adapted sparse grids are applied to effectively explore the parameter domain even if the number of parameters is modest or high. Adaptivity is also employed to build rational approximations that efficiently capture the frequency dependence of the problem. These two features enable our method to build surrogate models that achieve a user-prescribed approximation accuracy, without wasting too many resources in ``oversampling'' the frequency and parameter domains. Thanks to its non-intrusiveness, our proposed method, as opposed to projection-based techniques for model order reduction, can be applied regardless of the complexity of the underlying physical model. Notably, our algorithm for adaptive sampling can be used even when prior knowledge of the problem structure is not available. To showcase the effectiveness of our approach, we apply it in the study of an aerodynamic bearing. Our method allows us to build surrogate models that adequately identify the bearing's behavior with respect to both design and operational parameters, while still achieving significant speedups.