  Phase reduction theory has been applied to many systems with limit cycles;
however, it has limited applications in incompressible fluid systems. This is
because the calculation of the phase sensitivity function, one of the
fundamental functions in phase reduction theory, has a high computational cost
for systems with a large degree of freedom. Furthermore, incompressible fluid
systems have an implicit expression of the Jacobian. To address these issues,
we propose a new algorithm to numerically calculate the phase sensitivity
function. This algorithm does not require the explicit form of the Jacobian
along the limit cycle, and the computational time is significantly reduced,
compared with known methods. Along with the description of the method and
characteristics, two applications of the method are demonstrated. One
application is the traveling pulse in the FitzHugh Nagumo equation in a
periodic domain and the other is the K\'arm\'an's vortex street. The response
to the perturbation added to the K\'arm\'an's vortex street is discussed in
terms of both phase reduction theory and fluid mechanics.
