Shape memory alloys (SMAs) are a branch of smart materials, which have unique properties, including two important features: pseudoelasticity and shape memory effect. Shape memory effect is the recovery of permanent strains due to heating. Pseudoelasticity is the ability of material to recover large strains after the unloading process. Magnitude of these recoverable strains is up to 10%. These two important properties are due to alloy structure and two solid state phases: austenite and martensite. Austenite phase has high symmetry and high temperature and has more energy than martensite phase, which has low symmetry and low temperature. Martensitic transformation between these two phases occurs due to heatting or stress. SMA's have also variety of properties such as high power to weight ratio, energy damping and biocompatibility and are used in various industrial applications. These wide applications lead the necessity of modeling and simulation of SMA's. Different constitutive models have been presented for these alloys and different numerical approaches for the implementation of these models have been used. Most of the existing models are 1D and their solutions are only applicable for SMA wires lonely or for particular simple types of structures consisting SMA wires. Since smart structures consisting SMAs are various and have general geometries, 3D modeling and solutions is a necessity. In the last decades some 3D models have been provided. Most of 3D models are used for special geometries and load conditions. These special case studies are no longer effective in general cases. Coupled thermomechanical behavior of these materials is a challenge in modeling. This behavior occurs in high strain rates and dynamic loadings. In this situation the isothermal assumption is not valid and is needed to take account the temperature changes. In present research, an approach for modeling and implementing fully coupled thermomechanical behavior under multiaxial loadings is provided. This model includes both shape memory effect and pseudoelasticity and is able to simulates any thermal and multiaxial mechanical loading path. Also present model takes account the thermomechanical coupling and rate dependency of SMAs. The model implementation is carried out by ABAQUS commercial finite element package. Constitutive and energy relations are introduced into ABAQUS via a user material subroutine (UMAT). The model can investigate different problems with various boundary conditions, geometries and loadings. In order to validate the model and UMAT, model results are compared with results obtained from experimental, analytical and numerical efforts reported in the references. The comparison shows a good agreement in model results. The effect of different loading rates is investigated. Also the effects of thermal boundary conditions on temperature gradient are observed. Then, the effect of 3D geometry in a bar is investigated, and the results show the temperature gradient in cross section of bar and thus a different stress-strain curve in bar center and its surface is observed. As a case study helical springs are investigated and the model shows a good agreement in quasi-static loadings. In high loading rate, the behavior of helical spring is different from that of quasi-static condition. The results show that in loadings with high strain rate, latent heat is a key parameter. In this research the effect of latent heat is take into account explicitly. Keywords: Shape Memory alloy (SMA), Multiaxial loading, Fully coupled thermomechanical behavior, Latent heat, Finite element, ABAQUS, UMAT.