Uncertainties in the modelling of engineering systems arise, among others, from imprecisely known or varying material properties, manufacturing deviations or production tolerances, as well as different operating conditions. Accounting for them in systems with viscoelastic elements, treated as components with enhanced damping capacity, is crucial, because even small deviations may significantly alter the dynamic response of the structure. This is extremely important from the viewpoint of designing structures with controlled energy dissipation. The analysed systems include viscoelastic dampers and viscoelastic layers [1] described by classical and fractional rheological models. Uncertainties in design parameters are described in a non-probabilistic framework using interval analysis [2]. This approach is particularly useful when only limited identification data are available, which prevents a reliable determination of probabilistic distributions. To analyse uncertainty propagation, Kriging metamodeling is employed. The response of the analysed system, expressed in terms of dynamic characteristics, is treated as a function of multiple variables. Its approximation is determined based on cuts of the parameter space passing through a reference point, in accordance with the idea of High-Dimensional Model Representation (HDMR) [3]. Subsequently, Kriging metamodels are constructed; optimisation of these surrogates provides the lower and upper bounds of the analysed characteristics and the corresponding values of the uncertain parameters. In order to reduce the influence of the choice of the reference point, an iterative procedure is applied, consisting in updating the reference point on the basis of preliminarily identified extrema and reconstructing the metamodels. The proposed approach can support the development of reliable computational models of structures with high damping capacity under limited information on design parameters.