The use of fiber reinforced polymer (FRP) composites in strengthening reinforced concrete beam-column subassemblies has been scrutinised both experimentally and numerically in recent years. While a multitude of numerical models are available, and many match the experimental results reasonably well, there are not many studies that have looked at the efficiency of different finite elements in a comparative way in order to clearly identify the best practice when it comes to modelling FRP for strengthening. The present study aims at investigating this within the context of FRP retrofitted reinforced concrete beam-column subassemblies. Two programs are used side by side; ANSYS and VecTor2. Results of the finite element modeling using these two programs are compared with a recent experimental study. Different failure and yield criteria along with different element types are implemented and a useful technique, which can reduce the number of elements considerably, is successfully employed for modeling planar structures subjected to in-plane loading in ANSYS. Comparison of the results shows that there is good agreement between ANSYS and VecTor2 results in monotonic loading. However, unlike VecTor2 program, implicit version of ANSYS program is not able to properly model the cyclic behavior of the modeled subassemblies. The paper will be useful to those who wish to study FRP strengthening applications numerically as it provides an insight into the choice of the elements and the methods of modeling to achieve desired accuracy and numerical stability, a matter not so clearly explored in the past in any of the published literature.