A search for the Higgs boson decaying into a pair of charm quarks is presented. The analysis uses proton-proton collisions to target the production of a Higgs boson in association with a leptonically decaying $W$ or $Z$ boson. The dataset delivered by the LHC at a centre-of-mass energy of $\sqrt{s}=13$ TeV and recorded by the ATLAS detector corresponds to an integrated luminosity of 139 fb$^{-1}$. Flavour-tagging algorithms are used to identify jets originating from the hadronisation of charm quarks. The analysis method is validated with the simultaneous measurement of $WW,WZ$ and $ZZ$ production, with observed (expected) significances of $2.6$ ($2.2$) standard deviations above the background-only prediction for the $(W/Z)Z(\to c\bar{c})$ process and $3.8$ ($4.6$) standard deviations for the $(W/Z)W(\to cq)$ process. The $(W/Z)H(\to c \bar{c})$ search yields an observed (expected) upper limit of $26$ ($31$) times the predicted Standard Model cross-section times branching fraction for a Higgs boson with a mass of 125 GeV, corresponding to an observed (expected) constraint on the charm Yukawa coupling modifier $|\kappa_c| < 8.5 (12.4)$, at the 95% confidence level. A combination with the ATLAS $(W/Z)H, H\to b\bar{b}$ analysis is performed, allowing the ratio $\kappa_c / \kappa_b$ to be constrained to less than 4.5 at the 95% confidence level, smaller than the ratio of the $b$- and $c$-quark masses, and therefore determines the Higgs-charm coupling to be weaker than the Higgs-bottom coupling at the 95% confidence level.