  We develop a formalism for calculating forces on the nuclei within the
linear-scaling stochastic density functional theory (sDFT) in a nonorthogonal
atom-centered basis-set representation (Fabian et al. WIREs Comput Mol Sci.
2019;e1412. https://doi.org/10.1002/wcms.1412) and apply it to Tryptophan
Zipper 2 (Trp-zip2) peptide solvated in water. We use an embedded-fragment
approach to reduce the statistical errors (fluctuation and systematic bias),
where the entire peptide is the main fragment and the remaining 425 water
molecules are grouped into small fragments. We analyze the magnitude of the
statistical errors in the forces and find that the systematic bias is of the
order of $0.065\,eV/\r{A}$ ($\sim1.2\times10^{-3}E_{h}/a_{0}$) when 120
stochastic orbitals are used, independently of systems size. This magnitude of
bias is sufficiently small to ensure that the bond lengths estimated by
stochastic DFT (within a Langevin molecular dynamics simulation) will deviate
by less than 1% from those predicted by a deterministic calculation.
