The evaluation of contaminated sites with respect to their potential for natural attenuation to relevant orders of magnitude in concentration and reasonable time scales requires estimating future plume behaviour. In this regard, model simulations are a valuable tool, e. g. through studies of predictive scenarios, but they have to fulfil demanding criteria that common software codes often cannot meet. We first present an overview of work completed within the framework of the joint research project BayFoNA in the domain of modelling. The herein described universal software tool is capable of reflecting the relevant processes of contaminant propagation and transformation, including, for example carrier facilitation or general temperature dependent decay reactions. Due to site heterogeneity, highly localized reactions, as well as aspects including nonlinearity or varying time scales, the quality of the numerical methods can have a decisive influence on the assessment of the natural attenuation potential. Therefore non standard locally mass conservative finite element discretizations were applied, employing functions of higher order to reduce numerical dispersion (which, in some cases, may lead to an overestimation of microbial degradation), or a globally implicit solution algorithm to avoid splitting errors and convergence problems. The simulation tool also allows the identification of unknown parameters by well defined experiments. Applications in the joint research project demonstrate the implementation of the techniques in the case of BTEX contamination, and helped clarify issues related to the migration of arsenic.