A measurement of the four-lepton invariant mass spectrum is made with the ATLAS detector, using an integrated luminosity of 36.1 fb$^{-1}$ of proton-proton collisions at $\sqrt{s}$ = 13 TeV delivered by the Large Hadron Collider. The differential cross-section is measured for events containing two same-flavour opposite-sign lepton pairs. It exhibits a rich structure, with different mass regions dominated in the Standard Model by single $Z$ boson production, Higgs boson production, and $Z$ boson pair production, and non-negligible interference effects at high invariant masses. The measurement is compared with state-of-the-art Standard Model calculations, which are found to be consistent with the data. These calculations are used to interpret the data in terms of $gg\rightarrow ZZ \rightarrow 4\ell$ and $Z \rightarrow 4\ell$ subprocesses, and to place constraints on a possible contribution from physics beyond the Standard Model.

A search for the dimuon decay of the Higgs boson was performed using data corresponding to an integrated luminosity of 36.1 fb$^{-1}$ collected with the ATLAS detector in $pp$ collisions at $\sqrt{s}=13$ TeV at the Large Hadron Collider. No significant excess is observed above the expected background. The observed (expected) upper limit on the cross section times branching ratio is 3.0 (3.1) times the Standard Model prediction at the 95% confidence level for a Higgs boson mass of 125 GeV. When combined with the $pp$ collision data at $\sqrt{s}=7$ TeV and $\sqrt{s}=8$ TeV, the observed (expected) upper limit is 2.8 (2.9) times the Standard Model prediction.

This Letter presents direct searches for lepton flavour violation in Higgs boson decays, $H\rightarrow e\tau$ and $H\rightarrow\mu\tau$, performed with the ATLAS detector at the LHC. The searches are based on a data sample of proton-proton collisions at a centre-of-mass energy $\sqrt{s} = 13$ TeV, corresponding to an integrated luminosity of $36.1\,\mathrm{fb}^{-1}$. No significant excess is observed above the expected background from Standard Model processes. The observed (median expected) 95 % confidence-level upper limits on the lepton-flavour-violating branching ratios are $0.47\%$ ($0.34^{+0.13}_{-0.10}\,\%$) and $0.28\%$ ($0.37^{+0.14}_{-0.10}\,\%$) for $H\to e\tau$ and $H\to\mu\tau$, respectively.

A search for heavy neutral Higgs bosons is performed using the LHC Run 2 data, corresponding to an integrated luminosity of 139 fb$^{-1}$ of proton-proton collisions at $\sqrt{s}=13$ TeV recorded with the ATLAS detector. The search for heavy resonances is performed over the mass range 0.2-2.5 TeV for the $\tau^+\tau^-$ decay with at least one $\tau$-lepton decaying into final states with hadrons. The data are in good agreement with the background prediction of the Standard Model. In the $M_{h}^{125}$ scenario of the Minimal Supersymmetric Standard Model, values of $\tan\beta>8$ and $\tan\beta>21$ are excluded at the 95% confidence level for neutral Higgs boson masses of 1.0 TeV and 1.5 TeV, respectively, where $\tan\beta$ is the ratio of the vacuum expectation values of the two Higgs doublets.

A search for the $Z\gamma$ decay of the Higgs boson, with $Z$ boson decays into pairs of electrons or muons is presented. The analysis uses proton$-$proton collision data at $\sqrt{s}$ = 13 TeV corresponding to an integrated luminosity of 139 fb$^{-1}$ recorded by the ATLAS detector at the Large Hadron Collider. The observed data are consistent with the expected background with a $p$-value of 1.3%. An upper limit at 95% confidence level on the production cross-section times the branching ratio for $pp\to H\to Z\gamma$ is set at 3.6 times the Standard Model prediction while 2.6 times is expected in the presence of the Standard Model Higgs boson. The best-fit value for the signal yield normalised to the Standard Model prediction is $2.0^{+1.0}_{-0.9}$ where the statistical component of the uncertainty is dominant.

Measurements of the Standard Model Higgs boson decaying into a $b\bar{b}$ pair and produced in association with a $W$ or $Z$ boson decaying into leptons, using proton-proton collision data collected between 2015 and 2018 by the ATLAS detector, are presented. The measurements use collisions produced by the Large Hadron Collider at a centre-of-mass energy of $\sqrt{s} = $13 TeV, corresponding to an integrated luminosity of 139 fb$^{-1}$. The production of a Higgs boson in association with a $W$ or $Z$ boson is established with observed (expected) significances of 4.0 (4.1) and 5.3 (5.1) standard deviations, respectively. Cross-sections of associated production of a Higgs boson decaying into bottom quark pairs with an electroweak gauge boson, $W$ or $Z$, decaying into leptons are measured as a function of the gauge boson transverse momentum in kinematic fiducial volumes. The cross-section measurements are all consistent with the Standard Model expectations, and the total uncertainties vary from 30% in the high gauge boson transverse momentum regions to 85% in the low regions. Limits are subsequently set on the parameters of an effective Lagrangian sensitive to modifications of the $WH$ and $ZH$ processes as well as the Higgs boson decay into $b\bar{b}$.

A search for the Higgs boson decaying into a photon and a pair of electrons or muons with an invariant mass $m_{\ell\ell} < 30$ GeV is presented. The analysis is performed using 139 fb$^{-1}$ of proton-proton collision data, produced by the LHC at a centre-of-mass energy of 13 TeV and collected by the ATLAS experiment. Evidence for the $H \rightarrow \ell \ell \gamma$ process is found with a significance of 3.2$\sigma$ over the background-only hypothesis, compared to an expected significance of 2.1$\sigma$. The best-fit value of the signal strength parameter, defined as the ratio of the observed signal yield to the one expected in the Standard Model, is $\mu = 1.5 \pm 0.5$. The Higgs boson production cross-section times the $H \rightarrow\ell\ell\gamma$ branching ratio for $m_{\ell\ell} <$ 30 GeV is determined to be 8.7 $^{+2.8}_{-2.7}$ fb.

A measurement of the production processes of the recently discovered Higgs boson is performed in the two-photon final state using 5.4 fb$^{-1}$ of proton-proton collisions data at $\sqrt{s}=7$ TeV and 20.3 fb$^{-1}$ at $\sqrt{s}=8$ TeV collected by the ATLAS detector at the Large Hadron Collider. The number of observed Higgs boson decays to diphotons divided by the corresponding Standard Model prediction, called the signal strength, is found to be $\mu = 1.17 \pm 0.27$ at the value of the Higgs boson mass measured by ATLAS, $m_{H}$ = 125.4 GeV. The analysis is optimized to measure the signal strengths for individual Higgs boson production processes at this value of $m_{H}$. They are found to be $\mu_{\mathrm{ggF}} = 1.32 \pm 0.38$, $\mu_{\mathrm{VBF}} = 0.8 \pm 0.7$, $\mu_{{WH}} = 1.0 \pm 1.6 $, $\mu_{{ZH}} = 0.1 ^{+3.7}_{-0.1} $, $\mu_{{t\bar{t}H}} = 1.6 ^{+2.7}_{-1.8} $, for Higgs boson production through gluon fusion, vector-boson fusion, and in association with a $W$ or $Z$ boson or a top-quark pair, respectively. Compared with the previously published ATLAS analysis, the results reported here also benefit from a new energy calibration procedure for photons and the subsequent reduction of the systematic uncertainty on the diphoton mass resolution. No significant deviations from the predictions of the Standard Model are found.

A search for pair production of the supersymmetric partners of the Higgs boson (higgsinos $\tilde{H}$) in gauge-mediated scenarios is reported. Each higgsino is assumed to decay to a Higgs boson and a gravitino. Two complementary analyses, targeting high- and low-mass signals, are performed to maximize sensitivity. The two analyses utilize LHC $pp$ collision data at a center-of-mass energy $\sqrt{s} = 13$ TeV, the former with an integrated luminosity of 36.1 fb$^{-1}$ and the latter with 24.3 fb$^{-1}$, collected with the ATLAS detector in 2015 and 2016. The search is performed in events containing missing transverse momentum and several energetic jets, at least three of which must be identified as $b$-quark jets. No significant excess is found above the predicted background. Limits on the cross-section are set as a function of the mass of the $\tilde{H}$ in simplified models assuming production via mass-degenerate higgsinos decaying to a Higgs boson and a gravitino. Higgsinos with masses between 130 and 230 GeV and between 290 and 880 GeV are excluded at the 95% confidence level. Interpretations of the limits in terms of the branching ratio of the higgsino to a $Z$ boson or a Higgs boson are also presented, and a 45% branching ratio to a Higgs boson is excluded for $m_{\tilde{H}} \approx 400$ GeV.

A low-background inclusive search for new physics in events with same-sign dileptons is presented. The search uses proton-proton collisions corresponding to 20.3 fb$^{-1}$ of integrated luminosity taken in 2012 at a centre-of-mass energy of 8 TeV with the ATLAS detector at the LHC. Pairs of isolated leptons with the same electric charge and large transverse momenta of the type $e^{\pm}e^{\pm}, e^{\pm}\mu^{\pm}$, and $\mu^{\pm}\mu^{\pm}$ are selected and their invariant mass distribution is examined. No excess of events above the expected level of Standard Model background is found. The results are used to set upper limits on the cross-sections for processes beyond the Standard Model. Limits are placed as a function of the dilepton invariant mass within a fiducial region corresponding to the signal event selection criteria. Exclusion limits are also derived for a specific model of doubly charged Higgs boson production.

A search for Higgs boson pair production $pp \to hh$ is performed with 19.5 fb$^{-1}$ of proton--proton collision data at $\sqrt{s}=$ 8 TeV, which were recorded by the ATLAS detector at the Large Hadron Collider in 2012. The decay products of each Higgs boson are reconstructed as a high-momentum $b\bar{b}$ system with either a pair of small-radius jets or a single large-radius jet, the latter exploiting jet substructure techniques and associated $b$-tagged track-jets. No evidence for resonant or non-resonant Higgs boson pair production is observed. The data are interpreted in the context of the Randall--Sundrum model with a warped extra dimension as well as the two-Higgs-doublet model. An upper limit on the cross-section for $pp \to G^{*}_{\mathrm{KK}} \to hh \to b\bar{b} b\bar{b}$ of 3.2 (2.3) fb is set for a Kaluza--Klein graviton $G^{*}_{\mathrm{KK}}$ mass of 1.0 (1.5) TeV, at the 95\% confidence level. The search for non-resonant Standard Model $hh$ production sets an observed 95\% confidence level upper limit on the production cross-section $\sigma(pp \to hh \to b\bar{b}b\bar{b})$ of 202 fb, compared to a SM prediction of $\sigma(pp \to hh \to b\bar{b}b\bar{b}) = 3.6 \pm 0.5$ fb.

Measurements of the production cross-sections of the Standard Model (SM) Higgs boson ($H$) decaying into a pair of $\tau$-leptons are presented. The measurements use data collected with the ATLAS detector from $pp$ collisions produced at the Large Hadron Collider at a centre-of-mass energy of $\sqrt{s}=13\,\text{TeV}$, corresponding to an integrated luminosity of $139\,\text{fb}^{-1}$. Leptonic ($\tau\to\ell\nu_{\ell}\nu_{\tau}$) and hadronic ($\tau\to\text{hadrons}~\nu_{\tau}$) decays of the $\tau$-lepton are considered. All measurements account for the branching ratio of $H\to\tau\tau$ and are performed with a requirement $|y_H|<2.5$, where $y_H$ is the true Higgs boson rapidity. The cross-section of the $pp\to H\to\tau\tau$ process is measured to be $2.94 \pm 0.21 \text{(stat)} ^{+\,0.37}_{-\,0.32} \text{(syst)}$ pb, in agreement with the SM prediction of $3.17\pm0.09~ \mbox{pb}$. Inclusive cross-sections are determined separately for the four dominant production modes: $2.65 \pm 0.41 \text{(stat)} ^{+\,0.91}_{-\,0.67} \text{(syst)}$ pb for gluon$-$gluon fusion, $0.197 \pm 0.028 \text{(stat)} ^{+\,0.032}_{-\,0.026} \text{(syst)}$ pb for vector-boson fusion, $0.115 \pm 0.058 \text{(stat)} ^{+\,0.042}_{-\,0.040} \text{(syst)}$ pb for vector-boson associated production, and $0.033 \pm 0.031 \text{(stat)} ^{+\,0.022}_{-\,0.017} \text{(syst)}$ pb for top-quark pair associated production. Measurements in exclusive regions of the phase space, using the simplified template cross-section framework, are also performed. All results are in agreement with the SM predictions.

A search is presented for the direct pair production of a chargino and a neutralino $pp\to\tilde{\chi}^\pm_1\tilde{\chi}^0_2$, where the chargino decays to the lightest neutralino and the $W$ boson, $\tilde{\chi}^\pm_1 \to \tilde{\chi}^0_1 (W^{\pm}\to\ell^{\pm}\nu)$, while the neutralino decays to the lightest neutralino and the 125 GeV Higgs boson, $\tilde{\chi}^0_2 \to \tilde{\chi}^0_1 (h\to bb/\gamma\gamma/\ell^{\pm}\nu qq)$. The final states considered for the search have large missing transverse momentum, an isolated electron or muon, and one of the following: either two jets identified as originating from bottom quarks, or two photons, or a second electron or muon with the same electric charge. The analysis is based on 20.3 fb$^{-1}$ of $\sqrt{s}=8$ TeV proton-proton collision data delivered by the Large Hadron Collider and recorded with the ATLAS detector. Observations are consistent with the Standard Model expectations, and limits are set in the context of a simplified supersymmetric model.

A search for an invisibly decaying Higgs boson or dark matter candidates produced in association with a leptonically decaying $Z$ boson in proton--proton collisions at $\sqrt{s} =$ 13 TeV is presented. This search uses 36.1 fb$^{-1}$ of data collected by the ATLAS experiment at the Large Hadron Collider. No significant deviation from the expectation of the Standard Model backgrounds is observed. Assuming the Standard Model $ZH$ production cross-section, an observed (expected) upper limit of 67% (39%) at the 95% confidence level is set on the branching ratio of invisible decays of the Higgs boson with mass $m_H = $ 125 GeV. The corresponding limits on the production cross-section of the $ZH$ process with the invisible Higgs boson decays are also presented. Furthermore, exclusion limits on the dark matter candidate and mediator masses are reported in the framework of simplified dark matter models.

This paper reports searches for heavy resonances decaying into $ZZ$ or $ZW$ using data from proton--proton collisions at a centre-of-mass energy of $\sqrt{s}=13$ TeV. The data, corresponding to an integrated luminosity of 36.1 fb$^{-1}$, were recorded with the ATLAS detector in 2015 and 2016 at the Large Hadron Collider. The searches are performed in final states in which one $Z$ boson decays into either a pair of light charged leptons (electrons and muons) or a pair of neutrinos, and the associated $W$ boson or the other $Z$ boson decays hadronically. No evidence of the production of heavy resonances is observed. Upper bounds on the production cross sections of heavy resonances times their decay branching ratios to $ZZ$ or $ZW$ are derived in the mass range 300--5000 GeV within the context of Standard Model extensions with additional Higgs bosons, a heavy vector triplet or warped extra dimensions. Production through gluon--gluon fusion, Drell--Yan or vector-boson fusion are considered, depending on the assumed model.

A search for Higgs boson pair production in the $b\bar{b}b\bar{b}$ final state is carried out with up to 36.1 $\mathrm{fb}^{-1}$ of LHC proton--proton collision data collected at $\sqrt{s}$ = 13 TeV with the ATLAS detector in 2015 and 2016. Three benchmark signals are studied: a spin-2 graviton decaying into a Higgs boson pair, a scalar resonance decaying into a Higgs boson pair, and Standard Model non-resonant Higgs boson pair production. Two analyses are carried out, each implementing a particular technique for the event reconstruction that targets Higgs bosons reconstructed as pairs of jets or single boosted jets. The resonance mass range covered is 260--3000 GeV. The analyses are statistically combined and upper limits on the production cross section of Higgs boson pairs times branching ratio to $b\bar{b}b\bar{b}$ are set in each model. No significant excess is observed; the largest deviation of data over prediction is found at a mass of 280 GeV, corresponding to 2.3 standard deviations globally. The observed 95% confidence level upper limit on the non-resonant production is 13 times the Standard Model prediction.

This paper describes a search for beyond the Standard Model decays of the Higgs boson into a pair of new spin-0 particles subsequently decaying into $b$-quark pairs, $H \rightarrow aa \rightarrow (b\bar{b})(b\bar{b})$, using proton-proton collision data collected by the ATLAS detector at the Large Hadron Collider at center-of-mass energy $\sqrt{s}=13$ TeV. This search focuses on the regime where the decay products are collimated and in the range $15 \leq m_a \leq 30$ GeV and is complementary to a previous search in the same final state targeting the regime where the decay products are well separated and in the range $20 \leq m_a \leq 60$ GeV. A novel strategy for the identification of the $a \rightarrow b\bar{b}$ decays is deployed to enhance the efficiency for topologies with small separation angles. The search is performed with 36 fb$^{-1}$ of integrated luminosity collected in 2015 and 2016 and sets upper limits on the production cross-section of $H \rightarrow aa \rightarrow (b\bar{b})(b\bar{b})$, where the Higgs boson is produced in association with a $Z$ boson.

A search for heavy neutral Higgs bosons and $Z^{\prime}$ bosons is performed using a data sample corresponding to an integrated luminosity of 36.1 fb$^{-1}$ from proton-proton collisions at $\sqrt{s}$ = 13 TeV recorded by the ATLAS detector at the LHC during 2015 and 2016. The heavy resonance is assumed to decay to $\tau^+\tau^-$ with at least one tau lepton decaying to final states with hadrons and a neutrino. The search is performed in the mass range of 0.2-2.25 TeV for Higgs bosons and 0.2-4.0 TeV for $Z^{\prime}$ bosons. The data are in good agreement with the background predicted by the Standard Model. The results are interpreted in benchmark scenarios. In the context of the hMSSM scenario, the data exclude $\tan\beta > 1.0$ for $m_A$ = 0.25 TeV and $\tan\beta > 42$ for $m_A$ = 1.5 TeV at the 95% confidence level. For the Sequential Standard Model, $Z^{\prime}_\mathrm{SSM}$ with $m_{Z^{\prime}} < 2.42$ TeV is excluded at 95% confidence level, while $Z^{\prime}_\mathrm{NU}$ with $m_{Z^{\prime}} < 2.25$ TeV is excluded for the non-universal $G(221)$ model that exhibits enhanced couplings to third-generation fermions.

A search for heavy resonances decaying into a pair of $Z$ bosons leading to $\ell^+\ell^-\ell^+\ell^-$ and $\ell^+\ell^-\nu\bar\nu$ final states, where $\ell$ stands for either an electron or a muon, is presented. The search uses proton proton collision data at a centre-of-mass energy of 13 TeV corresponding to an integrated luminosity of 36.1 fb$^{-1}$ collected with the ATLAS detector during 2015 and 2016 at the Large Hadron Collider. Different mass ranges for the hypothetical resonances are considered, depending on the final state and model. The different ranges span between 200 GeV and 2000 GeV. The results are interpreted as upper limits on the production cross section of a spin 0 or spin 2 resonance. The upper limits for the spin 0 resonance are translated to exclusion contours in the context of Type I and Type II two-Higgs-doublet models, while those for the spin 2 resonance are used to constrain the Randall Sundrum model with an extra dimension giving rise to spin 2 graviton excitations.

The associated production of a Higgs boson with a $W$ or $Z$ boson decaying into leptons and where the Higgs boson decays to a $b\bar{b}$ pair is measured in the high vector-boson transverse momentum regime, above 250 GeV, with the ATLAS detector. The analysed data, corresponding to an integrated luminosity of 139 fb$^{-1}$, were collected in proton-proton collisions at the Large Hadron Collider between 2015 and 2018 at a centre-of-mass energy of $\sqrt{s} = 13$ TeV. The measured signal strength, defined as the ratio of the measured signal yield to that predicted by the Standard Model, is $0.72 ^{+0.39}_{-0.36}$ corresponding to an observed (expected) significance of 2.1 (2.7) standard deviations. Cross-sections of associated production of a Higgs boson decaying into $b$ quark pairs with a $W$ or $Z$ gauge boson, decaying into leptons, are measured in two exclusive vector boson transverse momentum regions, 250-400 GeV and above 400 GeV, and interpreted as constraints on anomalous couplings in the framework of a Standard Model effective field theory.

A search for heavy resonances decaying into a pair of $Z$ bosons leading to $\ell^+\ell^-\ell'^+\ell'^-$ and $\ell^+\ell^-\nu\bar\nu$ final states, where $\ell$ stands for either an electron or a muon, is presented. The search uses proton-proton collision data at a centre-of-mass energy of 13 TeV collected from 2015 to 2018 that corresponds to the full integrated luminosity of 139 fb$^{-1}$ recorded by the ATLAS detector during Run 2 of the Large Hadron Collider. Different mass ranges spanning 200 GeV to 2000 GeV for the hypothetical resonances are considered, depending on the final state and model. In the absence of a significant observed excess, the results are interpreted as upper limits on the production cross section of a spin-0 or spin-2 resonance. The upper limits for the spin-0 resonance are translated to exclusion contours in the context of Type-I and Type-II two-Higgs-doublet models, and the limits for the spin-2 resonance are used to constrain the Randall--Sundrum model with an extra dimension giving rise to spin-2 graviton excitations.

A study of the anomalous couplings of the Higgs boson to vector bosons, including $CP$-violation effects, has been conducted using its production and decay in the WW channel. This analysis is performed on proton-proton collision data collected with the CMS detector at the CERN LHC during 2016-2018 at a center-of-mass energy of 13 TeV, and corresponds to an integrated luminosity of 138 fb$^{-1}$. The different-flavor dilepton (e$\mu$) final state is analyzed, with dedicated categories targeting gluon fusion, electroweak vector boson fusion, and associated production with a W or Z boson. Kinematic information from associated jets is combined using matrix element techniques to increase the sensitivity to anomalous effects at the production vertex. A simultaneous measurement of four Higgs boson couplings to electroweak vector bosons is performed in the framework of a standard model effective field theory. All measurements are consistent with the expectations for the standard model Higgs boson and constraints are set on the fractional contribution of the anomalous couplings to the Higgs boson production cross section.

The first measurement of the CP structure of the Yukawa coupling between the Higgs boson and $\tau$ leptons is presented. The measurement is based on data collected in proton-proton collisions at $\sqrt{s} =$ 13 TeV by the CMS detector at the LHC, corresponding to an integrated luminosity of 137 fb$^{-1}$. The analysis uses the angular correlation between the decay planes of $\tau$ leptons produced in Higgs boson decays. The effective mixing angle between CP-even and CP-odd $\tau$ Yukawa couplings is found to be $-$1 $\pm$ 19$^\circ$, compared to an expected value of 0 $\pm$ 21$^\circ$ at the 68.3% confidence level. The data disfavour the pure CP-odd scenario at 3.0 standard deviations. The results are compatible with predictions for the standard model Higgs boson.

The results of a search for Higgs boson pair (HH) production in the WW*WW*, WW*$\tau\tau$, and $\tau\tau\tau\tau$ decay modes are presented. The search uses 138 fb$^{-1}$ of proton-proton collision data recorded by the CMS experiment at the LHC at a center-of-mass energy of 13 TeV from 2016 to 2018. Analyzed events contain two, three, or four reconstructed leptons, including electrons, muons, and hadronically decaying tau leptons. No evidence for a signal is found in the data. Upper limits are set on the cross section for nonresonant HH production, as well as resonant production in which a new heavy particle decays to a pair of Higgs bosons. For nonresonant production, the observed (expected) upper limit on the cross section at 95% confidence level (CL) is 21.3 (19.4) times the standard model (SM) prediction. The observed (expected) ratio of the trilinear Higgs boson self-coupling to its value in the SM is constrained to be within the interval $-$6.9 to 11.1 ($-$6.9 to 11.7) at 95% CL, and limits are set on a variety of new-physics models using an effective field theory approach. The observed (expected) limits on the cross section for resonant HH production range from 0.18 to 0.90 (0.08 to 1.06) pb at 95% CL for new heavy-particle masses in the range 250-1000 GeV.

Measurements of the inclusive and differential fiducial cross sections of the Higgs boson are presented, using the $\tau$ lepton decay channel. The differential cross sections are measured as functions of the Higgs boson transverse momentum, jet multiplicity, and transverse momentum of the leading jet in the event if any. The analysis is performed using proton-proton data collected with the CMS detector at the LHC at a center-of-mass energy of 13 TeV and corresponding to an integrated luminosity of 138 fb$^{-1}$. These are the first differential measurements of the Higgs boson cross section in the final state of two $\tau$ leptons, and they constitute a significant improvement over measurements in other final states in events with a large jet multiplicity or with a Lorentz-boosted Higgs boson.

A search is presented for the decay of the 125 GeV Higgs boson (H) to a pair of new light pseudoscalar bosons (a), followed by the prompt decay of each a boson to a bottom quark-antiquark pair, H $\to$ aa $\to$$\mathrm{b\bar{b}b\bar{b}}$. The analysis is performed using a data sample of proton-proton collisions collected with the CMS detector at a center-of-mass energy of 13 TeV, corresponding to an integrated luminosity of 138 fb$^{-1}$. To reduce the background from standard model processes, the search requires the Higgs boson to be produced in association with a leptonically decaying W or Z boson. The analysis probes the production of new light bosons in a 15 $\lt$$m_\mathrm{a}$$\lt$ 60 GeV mass range. Assuming the standard model predictions for the Higgs boson production cross sections for pp $\to$ WH and ZH, model independent upper limits at 95% confidence level are derived for the branching fraction $\mathcal{B}$(H $\to$ aa $\to$ $\mathrm{b\bar{b}b\bar{b}}$). The combined WH and ZH observed upper limit on the branching fraction ranges from 1.10 for $m_\mathrm{a} =$ 20 GeV to 0.36 for $m_\mathrm{a} =$ 60 GeV, complementing other measurements in the $\mu\mu\tau\tau$, $\tau\tau\tau\tau$ and bb$\ell\ell$ ($\ell=$ $\mu$,$\tau$) channels.