  We propose two closely--related Lagrangian numerical methods for the
simulation of physical processes involving advection, reaction and diffusion.
The methods are intended to be used in settings where the flow is nearly
incompressible and the P\'eclet numbers are so high that resolving all the
scales of motion is unfeasible. This is commonplace in ocean flows. Our methods
consist in augmenting the method of characteristics, which is suitable for
advection--reaction problems, with couplings among nearby particles, producing
fluxes that mimic diffusion, or unresolved small-scale transport. The methods
conserve mass, obey the maximum principle, and allow to tune the strength of
the diffusive terms down to zero, while avoiding unwanted numerical dissipation
effects.
