  In this paper we present a procedure to compute reducible invariant tori and
their stable and unstable manifolds in stroboscopic Poincar\'e maps. The method
has two steps. In the first step we compute, by means of a quadratically
convergent scheme, the Fourier series of the torus, its Floquet transformation,
and its Floquet matrix. If the torus has stable and/or unstable directions, in
the second step we compute the Taylor-Fourier expansions of the corresponding
invariant manifolds up to a given order. The paper also discusses the case in
which the torus is highly unstable so that a multiple shooting strategy is
needed to compute the torus. If the order of the Taylor expansion of the
manifolds is fixed and N is the number of Fourier modes, the whole
computational effort (torus and manifolds) increases as O(N log N) q and the
memory required behaves as O(N). This makes the algorithm very suitable to
compute high-dimensional tori for which a huge number of Fourier modes are
needed. Besides, the algorithm has a very high degree of parallelism. The paper
includes examples where we compute invariant tori (of dimensions up to 5) of
quasi-periodically forced ODEs. The computations are run in a parallel computer
and its efficiency with respect to the number of processors is also discussed.
