Sophisticated modeling of the migration of sorbing radionuclides in compacted claystones is needed for supporting the safety analysis of deep geological repositories for radioactive waste, which requires robust modeling tools/codes. Here, a benchmark related to a long term laboratory scale diffusion experiment of cesium, a moderately sorbing radionuclide, through Opalinus clay is presented. The benchmark was performed with the following codes: CORE2DV5, Flotran, COMSOL Multiphysics, OpenGeoSys-GEM, MCOTAC and PHREEQC v.3. The migration setup was solved with two different conceptual models, i) a single-species model by using a look-up table for a cesium sorption isotherm and ii) a multi-species diffusion model including a complex mechanistic cesium sorption model. The calculations were performed for three different cesium boundary concentrations (10−3, 10−5, 10−7 mol / L) to investigate the models/codes capabilities taking into account the nonlinear sorption behavior of cesium. Generally, good agreement for both single- and multi-species benchmark concepts could be achieved, however, some discrepancies have been identified, especially near the boundaries, where code specific spatial (and time) discretization had to be improved to achieve better agreement at the expense of longer computation times. In addition, the benchmark exercise yielded useful information on code performance, setup options, input and output data management, and post processing options. Finally, the comparison of single-species and multi-species model concepts showed that the single-species approach yielded generally earlier breakthrough, because this approach accounts neither for cation exchange of Cs+ with K+ and Na+, nor K+ and Na+ diffusion in the pore water.