  Chemotaxis systems play a crucial role in modeling the dynamics of bacterial
and cellular behaviors, including propagation, aggregation, and pattern
formation, all under the influence of chemical signals. One notable
characteristic of these systems is their ability to simulate concentration
phenomena, where cell density undergoes rapid growth near specific
concentration points or along certain curves. Such growth can result in
singular, spiky structures and lead to finite-time blowups.
  Our investigation focuses on the dynamics of the Patlak-Keller-Segel
chemotaxis system and its two-species extensions. In the latter case, different
species may exhibit distinct chemotactic sensitivities, giving rise to very
different rates of cell density growth. Such a situation may be extremely
challenging for numerical methods as they may fail to accurately capture the
blowup of the slower-growing species mainly due to excessive numerical
dissipation.
  In this paper, we propose a hybrid finite-difference-particle (FDP) method,
in which a particle method is used to solve the chemotaxis equation(s), while
finite difference schemes are employed to solve the chemoattractant equation.
Thanks to the low-dissipation nature of the particle method, the proposed
hybrid scheme is particularly adept at capturing the blowup behaviors in both
one- and two-species cases. The proposed hybrid FDP methods are tested on a
series of challenging examples, and the obtained numerical results demonstrate
that our hybrid method can provide sharp resolution of the singular structures
even with a relatively small number of particles. Moreover, in the two-species
case, our method adeptly captures the blowing-up solution for the component
with lower chemotactic sensitivity, a feature not observed in other works.
