This paper presents a novel physics-constrained adaptive polynomial chaos expansion surrogate model for multi-output systems, which commonly arise in multiphysics approximation, structural reliability analysis, and uncertainty quantification. In various typical physical problems and representative engineering applications, physics-constrained polynomial chaos expansion (PC^2) has proven to be effective, compared to conventional data-driven approach, in reducing computational cost and data dependence while preserving physical interpretability through enforced constraints [1]. However, existing PC^2 frameworks are primarily developed for single output models. For systems with multiple outputs, the frequently repeated solution of the constrained optimization becomes computationally costly. To overcome this limitation, a multivariate sensitivity-adaptive construction of a shared PCE basis, based on the admissible neighbor expansion, is incorporated into the existing PC^2 framework. In this formulation, all output components are represented within a common polynomial basis, and the corresponding coefficient sets are determined simultaneously through a unified constrained least squares. An adaptive basis selection strategy driven by coefficient sensitivity is employed to identify polynomial terms with joint relevance across all outputs while maintaining the hierarchical structure of the polynomial order [2]. The proposed approach provides an efficient surrogate model for multi-output systems governed by physical constraints, significantly reducing the computational cost of repeated constrained optimizations. Representative numerical example and an engineering case study will be presented to illustrate the applicability and performance of the proposed approach. The case study is represented by a complex finite element model representing a tunnel including soil-structure interaction, where it is crucial to accurately estimate inner forces reflecting geotechnical uncertainties. [1] Novak L., Sharma H., Shields M., Physics-Informed Polynomial Chaos Expansions, Journal of Computational Physics, Vol. 506, 112926, 2024. [2] Loukrezis D., Diehl E., De Gersem H., Multivariate Sensitivity-adaptive Polynomial Chaos Expansion for High-dimensional Surrogate Modeling and Uncertainty Quantification, Applied Mathematical Modelling, Vol. 137, 115746, 2025.