In the last few years we have developed stellar model atmospheres which included effects of anomalous abundances and strong magnetic field. The full treatment of anomalous Zeeman splitting and polarized radiative transfer were introduced in the model atmosphere calculations for the first time. In this investigation we present results of modelling the atmosphere of one of the most extreme magnetic chemically peculiar stars, HD137509. This Bp SiCrFe star has the mean surface magnetic field modulus of about 29kG. We use the recent version of the line-by-line opacity sampling stellar model atmosphere code LLmodels, which incorporates the full treatment of Zeeman splitting of spectral lines, detailed polarized radiative transfer and arbitrary abundances. We compare model predictions with photometric and spectroscopic observations of the star, aiming to reach a self-consistency between the abundance pattern derived from high-resolution spectra and abundances used for model atmosphere calculation. Based on magnetic model atmospheres, we redetermined abundances and fundamental parameters of HD137509 using spectroscopic and photometric observations. This allowed us to obtain a better agreement between observed and theoretical parameters compared to non-magnetic models with individual or scaled-solar abundances. We confirm that the magnetic field effects should be taken into account in the stellar parameter determination and abundance analysis.