In this paper, a gradient-based iterative algorithm is proposed for finding the least-squares solutions of the following constrained generalized inverse eigenvalue problem: given $$X\in C^{n\times m}$$ X∈Cn×m , $$\Lambda =\mathrm{diag}(\lambda _1,\lambda _2,\ldots ,\lambda _m)\in C^{m\times m}$$ Λ=diag(λ1,λ2,…,λm)∈Cm×m , find $$A^*,B^*\in C^{n\times n}$$ A∗,B∗∈Cn×n , such that $$\Vert AX-BX\Lambda \Vert $$ ‖AX-BXΛ‖ is minimized, where $$A^*,B^*$$ A∗,B∗ are Hermitian–Hamiltonian except for a special submatrix. For any initial constrained matrices, a solution pair $$(A^*,B^*)$$ (A∗,B∗) can be obtained in finite iteration steps by this iterative algorithm in the absence of roundoff errors. The least-norm solution can be obtained by choosing a special kind of initial matrix pencil. In addition, the unique optimal approximation solution to a given matrix pencil in the solution set of the above problem can also be obtained. A numerical example is given to show the efficiency of the proposed algorithm.