Treatment strategy for multiple intracranial aneurysms is challenging, as in many cases the choice of the order in which to treat aneurysms is not based on high-quality evidence. We aimed to digitally simulate clinical scenarios of two different orders in which multiple intracranial aneurysms were treated and analyze changes in hemodynamics after first stage of treatment. We prospectively included patients with two intracranial aneurysms, with order of treatment difficult to determine based on clinical data alone. For each patient we prepared three models of arteries harboring aneurysms: with both aneurysms present and with one of them removed. Computational modeling of blood flow using fluid-structure interaction methodology was performed for each model. Hemodynamic parameters of aneurysm domes were compared between models with both aneurysm present, and models in with aneurysms were removed in changing order. In 25 included patients, the calculated hemodynamic parameters such as Time-Averaged Wall Shear Stress (0.46 ± 0.40 vs. 0.54 ± 0.44 Pa
p <
0.01) and surface vortex fraction (12.73% ± 7.92% vs. 14.26% ± 7.46%
p = 0.02) decreased after first stage of treatment, while Time-Averaged Wall Shear Stress Gradient (1.44 ± 0.41 vs. 1.34 ± 0.46 Pa
p = 0.04) and percentage of wall shear stress <
0.5 Pa (50.13% ± 33.01% vs. 44.08% ± 34.16%
p <
0.01) increased. Changes of wall shear stress in remaining aneurysm dome were independently correlated with dome-to-neck ratio of both removed and remaining aneurysms. Hemodynamics of untreated aneurysm worsens after first stage of treatment. Dome-to-neck ratio of both treated and untreated aneurysm was the strongest and independent predictor of that change.