The pseudospectra, a tool to study the behavior of systems associated with nonnormal matrices has been considered extremely useful in the last decades, for that reason there is a recent interest in its efficient computation. In the large scale setting different projection Krylov methods have been used to calculate the pseudospectra, rather than earlier approaches which require the application of SVD decomposition at each point of a grid or inverse Lanczos schemes that still requires previous Schur factorization, resulting prohibited for larger matrices. In this work, we investigate the practical applicability and the performance of a new scheme to approximate the pseudospectrum of large matrices on high performance architectures. We propose a domain decomposition technique to get simultaneously the smallest singular value of A - zI for each point z = xk + yki on an convenient grid. Numerical results, on test matrices from the literature, are encouraging and show a reduction in time of this scheme compared with its sequential counterpart. Additionally, scalability features are promising.