This article presents the HDG approximation as a solution to the unilateral contact problem, leveraging the regularization method and an iterative procedure for resolution. In our study, u represents the potential (displacement of the elastic body) and q represents the flux (the force exerted on the body). Our analysis establishes that the utilization of polynomials of degree k(k≥1)\\documentclass[12pt]{minimal} \\usepackage{amsmath} \\usepackage{wasysym} \\usepackage{amsfonts} \\usepackage{amssymb} \\usepackage{amsbsy} \\usepackage{mathrsfs} \\usepackage{upgreek} \\setlength{\\oddsidemargin}{-69pt} \\begin{document}$k (k \\ge 1)$\\end{document} leads to achieving an optimal convergence rate of order k+1\\documentclass[12pt]{minimal} \\usepackage{amsmath} \\usepackage{wasysym} \\usepackage{amsfonts} \\usepackage{amssymb} \\usepackage{amsbsy} \\usepackage{mathrsfs} \\usepackage{upgreek} \\setlength{\\oddsidemargin}{-69pt} \\begin{document}$k+1$\\end{document} in L2\\documentclass[12pt]{minimal} \\usepackage{amsmath} \\usepackage{wasysym} \\usepackage{amsfonts} \\usepackage{amssymb} \\usepackage{amsbsy} \\usepackage{mathrsfs} \\usepackage{upgreek} \\setlength{\\oddsidemargin}{-69pt} \\begin{document}$L^{2}$\\end{document}-norm for both u and q. Importantly, this optimal convergence is maintained irrespective of whether the domain is discretized through a structured or unstructured grid. The numerical results consistently align with the theoretical findings, underscoring the effectiveness and reliability of the proposed HDG approximation method for unilateral contact problems.