The understanding of fossil fields origin, topology, and stability is one of the corner stones of the stellar magnetism theory. On one hand, since they survive on secular time scales, they may modify the structure and the evolution of their host stars. On the other hand, they must have a complex stable structure since it has been demonstrated that the simplest purely poloidal or toroidal fields are unstable on dynamical time scales. In this context, the only stable stellar configurations found today are those resulting from numerical simulations by Braithwaite and collaborators who studied the evolution of an initial stochastic magnetic field, which relaxes with a selective decay of magnetic helicity and energy, on mixed stable configurations (poloidal and toroidal) that seem to be in equilibrium and then diffuse. In this talk, we report the semi-analytical investigation of such an equilibrium field in the axisymmetric case. We use variational methods, which describe selective decay of magnetic helicity and energy during MHD relaxation, and we identify a supplementary invariant due to the stable stratification of stellar radiation zones. This leads to states that generalize force-free Taylor's relaxation states studied in plasma laboratory experiments that become non force-free in the stellar case. Moreover, astrophysical applications are presented and the stability of obtained configurations is studied.