  This paper describes a search for dark photons ($\gamma_d$) in proton-proton
collisions at $\sqrt{s}=13$ TeV at the Large Hadron Collider (LHC). The dark
photons are searched for in the decay of Higgs bosons ($H \to \gamma\gamma_d$)
produced through the $ZH$ production mode. The transverse mass of the system,
made of the photon and the missing transverse momentum from the non-interacting
$\gamma_d$, presents a distinctive signature as it peaks near the Higgs boson
mass. The results presented use the total Run-2 integrated luminosity of 139
fb$^{-1}$ recorded by the ATLAS detector at the LHC. The dominant reducible
background processes are estimated using data-driven techniques. A Boosted
Decision Tree technique is adopted to enhance the sensitivity of the search. As
no excess is observed with respect to the Standard Model prediction, an
observed (expected) upper limit on the branching ratio BR$(H\to
\gamma\gamma_d)$ of 2.28$\%$ (2.82$^{+1.33}_{-0.84}\%$) is set at 95$\%$ CL for
massless $\gamma_d$. For massive dark photons up to 40 GeV, the observed
(expected) upper limits on BR$(H\to \gamma\gamma_d)$ at 95\% confidence level
is found within the [2.19,2.52]$\%$ ([2.71,3.11]$\%$) range.
