In general, decision support is one of the main purposes of model-based analysis of systems. Response surface methodology (RSM) is an optimization technique that has been applied frequently in practice, but few automated variants are currently available. In this paper, we propose the combination of RSM with numerical analysis methods to solve continuous time Markov chain models of class-based queueing systems (CBQ). We consider first- and second-order models in RSM to identify an optimal parameter configuration for CBQ as part of the differentiated service architecture. Among the many known numerical solution methods for large Markov chains, we consider a Gauss-Seidel solver with relaxation that relies on a hierarchical Kronecker representation as implemented in the APNN Toolbox. To effectively apply the proposed optimization methodology we determine a suitable configuration of RSM and compare the results with previous results for optimizing CBQ.