An analytic para-Hermitian matrix is commonly used to characterize wideband signal models in array signal processing. In this paper, an algorithm for computing the analytic solution of a linear system involving an analytic, point-wise positive-definite, para-Hermitian matrix is proposed. The algorithm is based on the iterative Krylov-subspace method and extends the well-known conjugate gradient (CG) algorithm to an infinite-dimensional Hilbert space that contains the desired analytic solution. The convergence of the resulting para-Hermitian CG (pHCG) algorithm is proved based on the linearity, self-adjointness, boundedness, and positive-definiteness properties of a convolution operator associated with the para-Hermitian matrix. The proposed algorithm converges iteratively to the solution within a predefined accuracy, without involving matrix inversion or decomposition. In addition to the pHCG, which works with sequences (or time-domain signals), an equivalent $\boldsymbol{z}$-domain algorithm, named $\boldsymbol{z}$-pHCG, is also proposed. The effectiveness of the proposed algorithm is illustrated through numerical examples, which also reveal the limitations of conventional approaches from the literature.
