The volumetric reduction characteristic of distinct topology optimization techniques is attractive to both academia and industry. This interest stems from financial reasons associated with the reduction in material consumption and environmental aspect, which involves the reduction of natural resources and pollutant emissions. This study addresses a topology optimization formulation for structural elements based on the coupling between the Level Set Method (LSM) and the isogeometric formulation of the Boundary Element Method (IGABEM). The former method describes the boundary evolution whereas the latter determines the mechanical fields. As in IGABEM, isogeometric descriptions have been used in LSM, which provide accuracy, robustness, solution generality, and direct communication between both methods. The optimization problem can be defined according to the augmented Lagrangian method. Thus, the coupling between the methods has been written by defining the normal velocity to the reference level set curve through the shape sensitivity of the augmented Lagrangian function. A heuristic topology modification addresses the deficiency of LSM in generating holes. This study clarifies the possibility of encountering numerical instability problems when using the studied coupling formulation, which might also occur if another penalty method had been adopted for driving the optimization. An effective treatment is proposed for this issue, showing its capability to reduce the dependence on the parameters of the heuristic topology modification criterion for the process success. The results verify its robustness by comparing them with those from other coupling formulations between LSM and IGABEM where results from the Solid Isotropic Material with Penalization (SIMP) are the reference.