A numerical model for the dynamics of thin toroidal magnetic flux tubes in accretion disks of young stars is presented. The motion of magnetic flux tubes along the disk rotation axis is modeled taking into account buoyancy and drag forces, as well as the pressure of the external large-scale magnetic field. Evolution in the disk plane is calculated considering internal magnetic tension, external gas and magnetic pressures, gravity of the central object and centrifugal force. The model includes the effect of heat exchange with the surrounding gas. The disk structure is determined using the the MHD model of accretion disks by Dudorov and Khaibrakhmanov. The numerical model is implemented in a modular code in the C++ programming language. The dynamics equations are solved using 4th-order Runge-Kutta and Gear methods with automatic step control. The developed algorithms are tested on problems with known analytical solutions. Calculations for the parameters of the protoplanetary disk TW Hya show that under various conditions, the dynamical and thermal evolution of magnetic flux tubes can manifest as effective removal of mass and magnetic energy into the disk atmosphere with subsequent generation of outflows, accumulation of magnetic energy near the disk surface with subsequent flare phenomena, and formation of an inhomogeneous magnetized corona of the disk.