With that kind of output dimensionality-I’m not sure if you’ll get anything in a reasonable amount of time. GP’s are powerful-but computationally expensive. The number of inversions you need to perform to produce the posterior scales very quickly in the number of dimensions-and posterior uncertainty around the n-th vector will similarly scale very quickly.
Here’s a great intro lecture: https://www.youtube.com/watch?v=ttgUJtVJthA
Example notebook where you can find more links: Multi-output Gaussian Processes: Coregionalization models using Hamadard product — PyMC example gallery