Phosphatidylserine (PS) is a negatively charged lipid type commonly found in eukaryotic membranes, where it interacts with proteins via nonspecific electrostatic interactions as well as via specific binding. Moreover, in the presence of calcium ions, PS lipids can induce membrane fusion and phase separation. Molecular details of these phenomena remain poorly understood, partly because accurate models to interpret the experimental data have not been available. Here we gather a set of previously published experimental NMR data of C-H bond order parameter magnitudes, |SCH|, for pure PS and mixed PS:PC (phosphatidylcholine) lipid bilayers and augment this data set by measuring the signs of SCH in the PS headgroup using S-DROSS solid-state NMR spectroscopy. The augmented data set is then used to assess the accuracy of the PS headgroup structures in, and the cation binding to, PS-containing membranes in the most commonly used classical molecular dynamics (MD) force fields including CHARMM36, Lipid17, MacRog, Slipids, GROMOS-CKP, Berger, and variants. We show large discrepancies between different force fields and that none of them reproduces the NMR data within experimental accuracy. However, the best MD models can detect the most essential differences between PC and PS headgroup structures. The cation binding affinity is not captured correctly by any of the PS force fields - an observation that is in line with our previous results for PC lipids. Moreover, the simulated response of the PS headgroup to bound ions can differ from experiments even qualitatively. The collected experimental data set and simulation results will pave the way for development of lipid force fields that correctly describe the biologically relevant negatively charged membranes and their interactions with ions. This work is part of the NMRlipids open collaboration project (nmrlipids.blogspot.fi).