Is there is a way to use rstan or rstanarm to implement a generalized least squares regression? If not is there any plan to add this in the future? I know many people using the ape package to run PGSL who would be interested in being able to take a Bayesian approach if available and accessible.
If the covariance matrix of the errors is known up to a scale factor, then you can just multiply the outcome and the predictors by a Cholesky factor before feeding them to rstanarm::stan_glm. If you meant something more like Feasible Generalized Least Squares where the covariance matrix of the errors is estimated, then no there is nothing like that currently in rstanarm but you could write the Stan program for it yourself and estimate it with rstan.
Being able to provide the covariance matrix directly would be useful in situations where it is rank deficient and the cholesky decomposition cannot be computed from the start.