We present a thermodynamic Gibbs free energy model for the finite element simulation of the coupled thermo-magneto-mechanical behavior of ferromagnetic shape memory alloys (FSMAs). Starting from a free energy model for the conventional shape memory effect, additional terms are included to take into account the magnetic anisotropy and the geometry-dependent magnetostatic energy. Different functions are considered for the strain dependence of the anisotropy energy in order to describe the experimentally found strong dependence of the anisotropy energy on the ratio of short and long crystallographic axis c/a. The resulting energy landscape is used to calculate the transition probabilities between three martensite variants and the austenite state under applied stress and external magnetic field. The magnetic shape memory effect is simulated for different loading conditions and sample geometries. We demonstrate the influence of the c/a dependence of the anisotropy energy as well as the influence of twinning strain and elastic modulus on the transition between martensite variants. The model calculations are compared with experimental results on Ni-Mn-Ga single crystals.