The paper is devoted to developing an algorithm for solving the inverse problem of finding three-dimensional distribution of magnetization in a rectangular parallelepiped. This problem is ill-posed and described by the integral Fredholm equation of the first kind. The problem is reduced to solving a system of linear algebraic equations. To solve this system, the biconjucate gradient method is used. The implemented optimized algorithm exploits the Toeplitz-block-Toeplitz structure of the system matrix to reduce the memory requirements and computing time. © 2023 American Institute of Physics Inc.. All rights reserved.