An approach is proposed to model hemodynamics in an arteriovenous malformation and its vascular environment during neurosurgical embolization. This approach is based on a combination of the filtration flow model for blood and embolic agent in the malformation and the hydraulic approximation for the vessels surrounding the malformation. The model is described mathematically by a system of integrodifferential hyperbolic equations. The parameters and functions included in the model are determined using clinical data from real patients. Based on this model, the problem of optimal control of multistage embolization was formulated and studied numerically in a special class of controls. Optimal embolization regimes were found for which there is good agreement between the calculated and clinical data. The proposed approach can be used to develop preoperative recommendations about the optimal surgical intervention tactics.