Introduction. Most distribution networks are unbalanced and therefore require a specific solution for load flow. There are many works on the subject in the literature, but they mainly focus on simple network configurations. Among the methods dedicated to this problem, one can refer to the load flow method based on the bus injection to branch current and branch current to bus voltage matrices. Problem. Although this method is regarded as simple and complete, its drawback is the difficulty in supporting the transformer model as well as its winding connection types. Nevertheless, the method requires the system per unit to derive the load flow solution. Goal. In the present paper, our concern is the implementation of distribution transformers in the modeling and calculation of load flow in unbalanced networks. Methodology. Unlike previous method, distribution transformer model is introduced in the topology matrices without simplifying assumptions. Particularly, topology matrices were modified to take into account all winding types of both primary and secondary sides of transformer that conserve the equivalent scheme of an ideal transformer in series with an impedance. In addition, the adopted transformer models overcome the singularity problem that can be encountered when switching from the primary to the secondary side of transformer and inversely. Practical value. The proposed approach was applied to various distribution networks such as IEEE 4-nodes, IEEE 13-nodes and IEEE 37-nodes. The obtained results validate the method and show its effectiveness.