An interaction model considering the local stress on a particle surface is developed based on a volume averaging technique. With a scope to apply to turbulence modulation caused by particles of diameter comparable to the Kolmogorov length scale, grid width for resolving vortical structures outside the boundary layer of the particles is set to be close to the particle diameter. The interaction force in the volume averaged momentum equation is modelled based on analytical solutions of two fundamental flows. For. the uniform. flow case, the nonlinear effect of the first -order term in a series expansion with respect to the particle Reynolds number is found to be essential for the anisotropic Eulerian distribution of the interaction force. For the shear flow case, the anisotropic distribution of the interaction force is also essential, and it is modelled based on the Stokes's solution. Considering that the length scale of the averaging volume is determined to be comparable to the grid width and the particle diameter, the residual stress, which originates from the volume averaging of the nonlinear term in the momentum equation, is also modelled based on an undisturbed linear shear flow. According to the test simulation using the interaction force and residual stress models of the fundamental flows, the anisotropic interaction force model essentially improves the representation of the flow field and the mechanical work, and the effect of the residual stress is found to be reasonably reproduced by the present model. (C) 2016 Elsevier Ltd. All rights reserved.