On the Mass Accretion Rates of Herbig Ae/Be Stars. Magnetospheric Accretion or Boundary Layer?

Understanding how young stars gain their masses through disk-to-star accretion is of paramount importance in astrophysics. It affects our knowledge about the early stellar evolution, the disk lifetime and dissipation processes, the way the planets form on the smallest scales, or the connection to macroscopic parameters characterizing star-forming regions on the largest ones, among others. In turn, mass accretion rate estimates depend on the accretion paradigm assumed. For low-mass T Tauri stars with strong magnetic fields there is consensus that magnetospheric accretion (MA) is the driving mechanism, but the transfer of mass in massive young stellar objects with weak or negligible magnetic fields probably occurs directly from the disk to the star through a hot boundary layer (BL). The intermediate-mass Herbig Ae/Be (HAeBe) stars bridge the gap between both previous regimes and are still optically visible during the pre-main sequence phase, thus constituting a unique opportunity to test a possible change of accretion mode from MA to BL. This review deals with our estimates of accretion rates in HAeBes, critically discussing the different accretion paradigms. It shows that although mounting evidence supports that MA may extend to late-type HAes but not to early-type HBes, there is not yet a consensus on the validity of this scenario versus the BL one. Based on MA and BL shock modeling, it is argued that the ultraviolet regime could significantly contribute in the future to discriminating between these competing accretion scenarios.


Introduction
The evolution of young stars comprises several stages from the initial collapse in a molecular cloud until they enter the main sequence (MS), when the central objects reach enough temperature to burn hydrogen [1]. Although massive young stellar objects (MYSOs, defined from the stellar mass that will lead to a final supernova collapse; M * > 10 M ) keep their optically thick envelopes during their fast early evolution, lower-mass stars show an optically visible pre-main sequence (PMS) phase. Optically visible PMS stars are commonly classified based on their stellar masses. The lower-mass objects are the T Tauri stars (TT; 0.1 < M * /M < 2), divided into "Classical" (CTT) or "Weak" (WTT) depending on whether there are signs of ongoing accretion onto the central star. A stellar mass ∼2.5 M and temperature ∼8500 K is the rough limit below which the sub-photospheric regions are still fully convective (although convective envelopes are still present up to ∼4 M ) and thus in principle capable of generating magnetic fields through the dynamo process (e.g., [2,3] and references therein). As we will see in this review, the presence or absence of stellar magnetic fields is fundamental to understand how young stars gain their masses through accretion. The higher mass counterparts of CTTs are the Herbig Ae/Be stars (HAeBe; 2 < M * /M < 10). Essentially, HAeBes are young stars (≤10 Myr) with spectral types A and B, showing Hα and other emission lines in their spectra, and a circumstellar disk associated with infrared excess on top of arXiv:2005.01745v1 [astro-ph.SR] 4 May 2020 the photospheric emission. The general properties of HAeBes have been discussed in detailed specific reviews [4,5], and the reader can consult the online slides, proceedings, and collections associated with more recent conferences devoted to these stars (e.g., [6]). 1,2 Several approaches that involve vastly different spatial scales are necessary to understand how material collapses from larger to smaller structures that will finally lead to the formation of the individual stars. However, in last term understanding star formation requires knowing how the circumstellar material actually accretes onto the stellar surface at scales 1 au. In turn, understanding stellar accretion may have implications on the way that "macroscopic" parameters like the star formation rate (SFR) are estimated [7,8], or even on the formation process of planets at the smallest scales [9]. It is presently accepted that basically all young stars are surrounded by circumstellar disks, and even MYSOs may accrete a non-negligible part of their final masses through these structures (see e.g., the review in [10]). Therefore, it is necessary to obtain accurate estimates of the disk-to-star mass accretion rate (Ṁ acc ) and thus to know how such an accretion proceeds across a large range of stellar masses. Indeed, derivingṀ acc values requires a formal scenario from which direct observations can be interpreted, and different paradigms can lead to different accretion rate estimates based on the same observational data. For CTTs there is consensus that accretion is magnetically driven according to the magnetospheric accretion scenario (MA [11][12][13]), while for more massive stars without magnetic fields accretion may proceed directly from the disk to the star through a boundary layer (BL [14]). In this respect, HAeBes represent a fundamental regime that bridges the gap between the accretion properties of low-mass CTTs and those of MYSOs. Moreover, early-type Herbig Be stars (HBes) are the most massive stars for which direct accretion signatures can still be observed, given that MYSOs embedded in their natal clouds are opaque to robust accretion tracers that emit in the optical and ultraviolet (UV).
This review focuses on our estimates of accretion rates in HAeBe stars, thus discussing the way that disk-to-star accretion may proceed in these sources. Given that much of our current understanding of this topic has been partially inspired by the better known TT stars, Section 2 starts with an historical overview about how accretion has been understood and measured first for these objects and then for the HAeBes. Section 3 critically discusses the viability of MA in HAeBes mainly focused on the required and the observed magnetic fields. The different ways to measure accretion rates and the corresponding accuracies based on MA are described in Section 4. In Section 5 the few accretion rate estimates based on the BL scenario that are available for the HAeBe regime are discussed in comparison with the MA measurements. Then it is argued in Section 6 that the UV regime may be critical to test the validity of both competing scenarios. Finally, Section 7 includes some concluding remarks.

A Brief Historical Perspective
Although HAeBes tend to be brighter than TTs, which facilitates the detailed study of some of their disks through high-spatial resolution techniques, TT stars are generally better understood than the intermediate-and high-mass star regimes. Apart from being the precursors of Solar-like stars similar to our own, the main reason is that TTs are comparatively easier to find. This results from the fact that the shape of the initial mass function favors the formation of low-mass objects, and because the PMS phase is shorter as the stellar mass increases. Thus, in many aspects-and concerning accretion in particular-our knowledge of HAeBes is at least partially guided by previous works on TTs. 1 http://www.eso.org/sci/meetings/2014/haebe2014.html. contradictory results (see the review in [46]). In particular, the accretion rates estimated from the IR excesses and the BL scenario by Hillebrand et al. (1992) [47] were in the range 10 −5 -10 −7 M yr −1 , which are too high compared to more recent estimates (see below, Section 5, and [48]). Sorelli et al. (1996) [49] then suggested that the redshifted absorptions observed in the NaI doublet lines of some HAes showing UXOr-like photo-polarimetric variability could be explained in the context of MA. A few years later, the first indications that MA could be a valid scenario for the HAes as a whole-but not for HBes-were based on spectropolarimetry. In fact, Vink et al. (2002) [50] found a difference between the Hα spectropolarimetric signatures of HAe and HBe stars, which they interpreted as a transition from magnetically driven to direct, disk-to-star accretion depending on the spectral type. Similar spectropolarimetric studies including more stars and spectral lines have confirmed this change of behavior, which has been related to different accretion modes in CTTs and HAe stars on the one hand and in HBes on the other [51][52][53]. Almost in parallel, Eisner et al. (2004) [54] also suggested that there may be a transition from MA in late-type HAeBes to disk accretion in early-type sources, this time based on different inner disk geometries as inferred from near-infrared (IR) interferometry. That same year, Calvet et al. (2004) [55] presented results based on MA shock modeling applied to a small sample of "intermediate-mass TT stars" (IMTTs) with properties in-between TTs and HAeBes. They showed not only that the near-UV excesses and line profiles of these sources are consistent with the MA paradigm, but also that the empirical correlation between the accretion luminosities from MA shock modeling and the luminosity of the Brγ emission line extends from CTTs to IMTTs. However, the first detailed MA line and shock modeling applied to reproduce the observations of a HAe star was in the seminal paper by Muzerolle et al. (2004) [56]. In this work several optical line profiles of the prototypical star UX Ori were reproduced from MA. Moreover, Muzerolle et al. (2004) [56] suggested that the near-UV excess flux observed in the Balmer region of the spectra (∼3000-4500 Å) of many HAeBes [57] could be explained from MA shock modeling in a similar way as for CTT stars, establishing a calibration relating the observed "Balmer excess" (∆D B , see the top left panel of Figure 2) in UX Ori withṀ acc .
The previous works were then extrapolated to infer initial estimates of accretion rates based on MA for relatively wide samples of HAeBes. In particular, the empirical correlation with the Brγ luminosity observed by Calvet et al. (2004) [55] was used by Garcia-Lopez et al. (2006) [58] to inferṀ acc values for dozens of HAes. Although Calvet et al. (2004) [55] derived the correlation from a sample of IMTTs, a very similar correlation was later found for the HAeBes (see below), for which the extrapolation by Garcia-Lopez et al. (2006) [58] proved accurate. Similarly, the ∆D B -Ṁ acc calibration by Muzerolle et al. (2004) [56] was also applied to estimate accretion rates of dozens of HAeBes with a wide range of stellar properties [59,60], although that calibration is only valid for stars with the same stellar parameters than UX Ori. As we will see next and in Section 6, the extrapolation of such a calibration to HAeBe stars with different stellar properties can lead to accretion rates systematically biased by more than an order of magnitude.
The first self-consistent estimates ofṀ acc based on MA for a wide sample of HAeBes were made by   [61] from the observed Balmer excesses of 38 northern stars. The photometric excess of each object was reproduced using the MA shock models of   [24], deriving individual mass accretion rates for the whole sample and demonstrating that the calibration ∆D B −Ṁ acc is strongly dependent on the specific stellar properties. In particular, a given ∆D B translates into significantly higher accretion rates as the stars are hotter, and especially as the stellar surface gravity decreases (i.e., for smaller M * /R * ratios). In addition, that paper showed that the widely used empirical calibrations between the accretion luminosity (L acc ) and the luminosity of emission lines in CTTs can be extended to the HAeBes, at least for the Hα, [OI]6300, and Brγ lines studied in that work (see also [62]). Interestingly,   [61] also reported that it is impossible to reproduce the strong Balmer excesses of a few HBe stars in their sample from MA shock modeling, suggesting that an alternative accretion mechanism may operate in these objects and constituting a strong support to the initial claims from spectropolarimetry and interferometry mentioned above. Later, Fairlamb et al. (2015) [63] applied a similar methodology to derive accretion rates from MA shock modeling and X-Shooter spectra of 91 southern HAeBes. Again, a significant fraction (>25%) of the HBe stars in that new sample could not be fitted from MA, reinforcing the view that the accretion physics could change for the stars with the earliest spectral types. The use of the X-Shooter spectra covering a wide wavelength range from the near-UV to the near-IR led Fairlamb et al. (2017) [64] to update previous empirical correlations from   [61] and to find new ones between L acc and the luminosity of many other emission lines, which are very similar to the corresponding correlations in CTTs [42]. The origin of these intriguing correlations in both CTTs and HAeBes was studied in   [65]. This work showed that indeed all lines from the near-UV to the near-IR can be used to infer accretion luminosities-even if some may not be physically related with accretion-because they reflect an underlying relation between L acc and the stellar luminosity, L * . In fact, the recent work by Wichittanakom et al. (2020) [66] provides an empirical calibration between L acc and L * for HAeBes that can also be used to derive rough estimates of averaged accretion rates. That the slope of the L acc -L * empirical calibration is shallower for HBes than for HAes and CTTs was interpreted by Wichittanakom et al. (2020) [66] as the signature of a different physical mechanism driving accretion; MA in HAes, and BL in HBes. Figure 1. Cut in the plane of the star perpendicular to the left side of an edge-on accreting disk where the dust (red) destruction radius is further from the star than the gas disk (blue). Three possible scenarios are shown corresponding to decreasing strengths of the stellar magnetic field from top to bottom. Gas is channeled through the field lines according to MA (top and middle panels, corresponding to decreasing sizes of the magnetosphere), and directly onto the star through a BL (indicated in cyan at the bottom panel) in the absence of a strong enough magnetic field.
Currently, Gaia has allowed the characterization of hundreds of known HAeBe stars [67] and accretion rate estimates from the empirical correlations with the Hα or the stellar luminosities are available for most of them [66,68,69]. The typical accretion rates of HAes from MA, ∼10 −7 M yr −1 , are an order of magnitude larger than for CTTs, and there is again great scatter and a scaling relation with the stellar mass.
It is noted that the fact that we can derive accretion rates of HAeBes based on MA does not mean that other paradigms cannot reproduce at least some of the observational properties of the same objects, and thus that they could also provide alternative accretion rates. Unfortunately, the BL view has received much less theoretical attention and only a few works provide accretion rate estimates of HAeBe stars based on this scenario [47,48]. The following section discusses the main arguments supporting and challenging MA as a main driving mechanism in HAeBes, and in Sections 5 and 6 the BL measurements and how they compare to MA estimates are discussed.

Is Magnetospheric Accretion Plausible for Herbig Ae/Be Stars?
So far we have seen that initial works suggested that the validity of MA could be extended at least to the late-type HAeBes, later supported by MA modeling that has led to direct measurements of accretion rates and to the extension of the empirical correlations as indirect probes of accretion. In addition, different lines of evidence support MA as the driving mechanism at least for HAe stars. For instance, multi-epoch analysis finds that the timescales and the extent of the Hα variability is similar for objects ranging in mass from 0.1 to 5 M [70] but smaller for HBe stars [71]. Statistical studies of wide samples of HAeBes in the optical and the near-IR show that the presence of redshifted and blueshifted self-absorptions in several emission lines are consistent with MA acting in HAes (although with small magnetospheres; see below) but not in HBes, for which the BL scenario seems more suitable [72,73]. Direct constraints of the Brγ and Hα line emitting regions from current spectro-interferometric facilities reveal that while those are in principle small enough to be consistent with the expected MA sizes in several HAeBes, many show more extended regions likely indicating an additional wind component (see e.g., [74] and references therein). Analogous spectro-interferometric studies devoted to fainter T Tauri stars are less frequent (see e.g., [75,76]), but in fact their emission line profiles are currently reproduced using hybrid models that consider magnetically driven accretion and winds (e.g., [33,34]). Similarly, the emission line profiles of HAeBes can also be reproduced either from MA modeling alone or combined with magnetically driven winds [56,61,74,[77][78][79][80][81]. Nonetheless, although many of these models are natural extensions from previous devoted to CTT stars, they lack a careful treatment of the larger projected rotational velocities that are characteristic of many HAeBes, which still constitutes a major limitation for the applicability of MA line modeling in these objects [28,74,81]. Other direct and indirect observational tests mentioned in the previous section that were important to support MA in CTTs can hardly be applied to HAeBes. In particular, hot accretion spots cannot be generally observed because their temperatures and the stellar ones are comparable for HAeBes [56]. Moreover, the potential influence of accretion on the stellar rotation is also difficult to measure given that non-accreting "weak HAeBes" analogous to the WTTs are not well identified (but see Section 4 and [61,82,83]).
On the other hand, a few works suggest that MA could not be valid for HAeBes. For instance, the analysis of multi-epoch spectroscopic observations of two HAeBe stars suggest that the absence of inverse P Cygni profiles and the stronger variability at longer timescales observed in the blueshifted absorptions of spectral lines are inconsistent with a scaled-up T Tauri MA scenario [84]. Similarly, that the fraction of redshifted and blueshifted absorption profiles in the HeI 10830 Å line of 5 "magnetic" and 59 "non-magnetic" HAeBe stars is similar may indicate that the stellar magnetic field does not play a role in the gas kinematics [85].
In fact, the primary requirement for MA to apply is the presence of a strong enough magnetic field capable of truncating the inner disk and channeling the material towards the stellar surface. Magnetic fields of the order of ∼1 kG are commonly detected in CTT stars, based on ≥ a dozen of such sources for which measurements have been carried out (see [3,22] and references therein). Those strong magnetic fields are enough to sustain MA in CTTs (e.g., [20][21][22]). On the opposite, the vast majority of HAeBes show magnetic fields ≤ hundreds of G or below detection limits, as inferred from wide samples including dozens of sources (e.g., [86][87][88][89][90][91]). The contrast between both previous results is probably behind the strongest arguments against MA working in HAeBes, and it certainly motivates the question on whether this scenario is plausible for the HAeBe regime or not (see e.g., the corresponding discussions in [63,72,81,88,92]).
Following Johns-Krull et al. (1999) [20], a lower limit of the magnetic field required to drive MA as a function of stellar and accretion parameters is given by: (see [93]). The stellar rotation period P * can be expressed as 2πR * sin i/vsin i, with vsin i the projected rotational velocity. Therefore, Equation (1) indicates that although the necessary B min increases with the stellar mass and accretion rate, larger for HAeBes than for CTTs, the major dependence is on the inverse of the stellar radius and rotational velocity. Given that both previous parameters tend to be substantially larger for more massive objects, the value of B min is typically smaller for HAeBes than for CTTs. In last term, this results from the fact that at distances below the co-rotation radius accretion dominates over the magnetic pressure driving winds, and such a radius tends to be comparatively smaller in HAeBes. Table 1 shows the compilation by Hubrig et al. (2015) [94] listing a representative sample of HAeBes with large-scale, organized magnetic fields averaged from detections made by different teams and instrumentation (B measured ). The stellar and accretion parameters necessary to derive the minimum magnetic fields from Equation (1) are also shown. The comparison between the last two columns shows that B measured ≥ B min in most cases (considering the errorbars of the measured values), as expected if MA is the driving mechanism. There are only four possible exceptions: BF Ori, HD 101412, HD 104237, and HD 190073. However, the use of slightly different stellar parameters and accretion rates could translate into smaller values of B min and make them consistent with the B measured values. For instance, the stellar parameters and accretion rate for BF Ori in   [61] lead to B min ∼50 G, smaller than the one measured and thus in potential agreement with MA too. Similarly, close binarity can introduce important uncertainties in the stellar parameters, accretion rates, and measured magnetic fields. In fact, HD 104237 is a spectroscopic binary with a HAe primary and a TT secondary and a recent work reports that the primary has a weak magnetic field ranging from ∼47 to 72 G [95]. HD 101412 shows the largest difference between B min and B measured , which is more than ∼ an order of magnitude smaller. Nonetheless, other measurements of the same star suggested that this could indeed be one of the HAeBes with the strongest magnetic field of several kG (see [96] and references therein), which combined with uncertainties in stellar parameters and accretion rates could be potentially consistent with the minimum magnetic field necessary to drive MA. Finally, HD 190073 is a massive HBe star and the most recent measurements in Järvinen et al. (2019) [91] confirm that its magnetic field is below ∼100 G, a factor >6 smaller than the value of B min estimated here and ruling out MA as the driving mechanism for this source.
Other analytical expressions making different assumptions of the coupling of the magnetic field with the inner disk yield constraints somewhat different than Equation (1) (see references in [20,21]). However, considering the uncertainties involved the typical stellar and accretion parameters of HAeBes lead to minimum magnetic fields that are comparable to the different measurements in the literature cited above; of the same order than the ones included in Table 1 for the representative sample of HAeBes with averaged detections from different teams. Concerning the HAeBes with non-detections, it must be emphasized that detection limits of several hundreds G are common [88]. While such limits are not problematic for TTs with kG magnetic fields, they are of the same order than the magnetic fields necessary to drive MA in most HAeBes. Moreover, the uncertainties involved and the magnetic field detection limits are generally larger for HAeBes than for TTs. This is because, first, magnetic field estimates are based on measurements of different photospheric lines that are less abundant in HAeBes than in TTs, reducing the statistical reliability. Perhaps more importantly, rotational line broadening constitutes an important limitation in measuring the strength and configuration of stellar magnetic fields. In particular, Villebrun et al. (2019) [3] recently showed that for vsin i > 80 km s −1 the detection limits could be large enough to prevent firm detections of large-scale magnetic fields, and such limits can be >1 kG for vsin i > 50 km s −1 for ordered, dipolar magnetic fields. In turn, HAeBes commonly show vsin i ≥ 100 km s −1 . Thus, it is unclear whether non-detections in HAeBes refer to an actual absence of a strong enough magnetic field to drive MA or not.  [85], except for BF Ori [61]; HD 58647 [97,98]; Z CMa [99,100]; HD 97048 [88,101]; HD 98922 (inclination from [84]); HD 100546 [102,103] and HD 104237 [104,105]. For HD 176386 an inclination of 50 • is assumed. Rotational periods are derived from the values for R max * , vsin i max and i min , except for V380 Ori [106].
In summary, although the absence of kG magnetic fields in HAeBes is commonly argued as a main flaw against MA working in these stars, such strong magnetic fields are necessary for lower-mass CTT stars but not for HAeBes, which require much smaller values to drive accretion magnetically. In fact, the minimum magnetic fields necessary to drive MA in HAeBes are generally consistent with current measurements and actual detection limits. In this respect, MA would remain as a plausible scenario also for HAeBes. Nonetheless, more work is still necessary concerning magnetic fields in HAeBes, in particular aiming at obtaining high signal-to-noise observations that reduce the detection limits and are capable of clarifying the somewhat contradictory results sometimes obtained by different teams and observational techniques (see e.g., [107]). The magnetic field issue actually remains as a fundamental observational and theoretical problem that affects our understanding of how intermediate-mass stars grow. In turn, if MA is the driving mechanism in HAeBes their typically smaller magnetic fields and larger rotational velocities naturally lead to more compact magnetospheres than for CTTs. The typical disk truncation radius for HAeBes would be ∼2.5 R * , either based on analytical formula depending on the magnetic field (e.g., [12]) or on the co-rotation radius (that sets a maximum radial distance from which gas can be magnetically driven onto the star; e.g., [13]). That radius is ∼ two times smaller than the one typically estimated for CTT stars (e.g., [24]). Eventually, the disk truncation radius could reach the stellar surface for small enough magnetic fields, when we would enter the BL regime (Figure 1).

Magnetospheric Accretion Measurements of HAeBe Stars
Although there is not yet a consensus about the physical validity of MA in HAeBe stars, this scenario is the one that has been more frequently used to derive accretion rates in these sources (but see Section 5). In this section, the MA-based approaches from which accretion rates can (or cannot) be derived for HAeBes and the corresponding accuracies are discussed (see also Section 5 in the review of [92] for a discussion on the lower limits).
As introduced in Section 2, accretion rates in HAeBes can be observationally estimated either from direct methods involving accretion shock or emission line modeling, or from indirect methods comprising empirical correlations with the previous, model-based accretion rates. Both the direct and indirect methods are exemplified in Figure 2. The direct methods are summarized in the following.

•
Accretion shock modeling ( Figure 2, top left): This is the most accurate method to directly infer a value ofṀ acc by reproducing the observed excess in the near-UV region of the spectrum from the contribution of the accretion shock. Such a contribution can be modeled as a blackbody for HAeBes [56,61], and its influence depends on free parameters relatively constrained by theory, like the inward flux of energy carried by the accretion funnels and the fraction of the stellar surface covered by the accretion shocks. The wavelength region of study in HAeBes has so far been centered on the Balmer jump, ∆D B , but it could be extended to shorter wavelengths (Section 6, and see [35] for TTs). The value of ∆D B , and thus the inferred value ofṀ acc , is itself distance-independent, but it could indirectly depend on the distance given that the measured excess depends on the stellar parameters assumed. Considering usual uncertainties and dependencies, the typical errorbar forṀ acc as estimated from this method is <0.5 dex. Details about the procedure, uncertainties involved, and measurements for dozens of northern and southern HAeBes can be consulted in the literature [61,63]. • Emission line modeling ( Figure 2, top right): According to MA several emission lines like the Balmer series, the sodium doublet, or HeI transitions are at least partially generated in the hot gas free-falling within the magnetic channels that connect the inner disk and the accretion shocks on the stellar surface. Radiative transfer is applied normally assuming a simple dipole geometry for the magnetic field. Apart from the stellar parameters (T * , M * , R * , and rotational velocity), MA line modeling depends on the disk inclination, the size of the magnetosphere-i.e., the disk truncation radius-and the gas temperature, on top ofṀ acc which is normally the parameter that one wants to determine. Even in the best case scenario when the stellar parameters and geometry are well constrained, so far we can only derive upper limits on the disk truncation radius based on spectro-interferometry, although more direct constraints can be inferred from spectro-polarimetry for a few stars (the potential of this latter technique to probe such small scales is discussed in [108,109]). In turn, rough estimates of the gas temperature are solely based on empirical constraints but theory is still lacking in this respect [28]. Given the number of free parameters, MA line modeling normally serves to estimateṀ acc within ∼ an order of magnitude accuracy, although for well-known sources and a careful modeling the uncertainty can be significantly reduced (see e.g., the recent work for a TT star in [110]). As noted in Section 3, MA line modeling including a complete treatment of the high rotational velocities that are typical in many HAeBes is still pending. Examples and details of line modeling applied to HAeBes, either considering MA alone or combined with magnetically driven winds, can be found in the literature [56,61,74,[77][78][79][80][81]. The indirect methods are discussed next: • Empirical correlations with the luminosity of emission lines (Figure 2, bottom left): The accretion luminosities of HAeBes inferred from the direct methods described above (L acc ∼ GM * Ṁacc /R * ) correlate with the luminosities of dozens of emission lines spreading from the near-UV to the near-IR, as previously found for CTTs (the extension of such correlations to specific lines at shorter and longer wavelengths in the far-UV and the mid-IR can be found e.g., in [111,112] respectively). Regardless of the physical origin of those lines all can be used to derive accretion rates by measuring the emission line luminosity and using the corresponding empirical expression with the form log (L acc /L ) = A( ± A) + B( ± B) × log (L line /L ) [65]. The errors for the slopes and intercepts are determined from the least-square fits to the L acc -L line data, and the final typical uncertainty forṀ acc (once transformed from L acc using the above mentioned formula and the stellar parameters) is ∼ 1 dex. This uncertainty could in principle be reduced by averaging the results obtained from different emission lines [41]. The most recent L acc -L line empirical expressions for HAeBes including dozens of emission lines can be found in Fairlamb et al. (2017) [64], and accretion rates inferred from this method (in particular, from the correlation with L Hα ) can also be found in the literature for hundreds of HAeBes [66,68].

•
Empirical correlations with the stellar luminosity ( Figure 2, bottom right): As mentioned above the fact that the luminosity of a given emission line correlates with L acc does not necessarily mean that there is an actual physical link between the origin of the line and the accretion process, and such a correlation naturally results from the underlying one between L acc and L * [65]. Because the scatter in the latter correlation is similar than for the L acc -L line correlations, ∼ ± 1 dex, accretion rates can be similarly derived from the L acc -L * empirical expression. This expression depends on the mass regime as described in Wichittanakom et al. (2020) [66] including HAeBes, and has been recently used to derive accretion rates for almost all HAeBes known [69].
Concerning the indirect methods, these do not generally serve to probe accretion variability except strong changes over relatively long timescales. The intrinsic scatter of ∼ ±1 dex that affects the L acc -L line correlations limits the range of accretion rate variations that could be analyzed. Still, L line variations homogeneously measured for a given star could be potentially linked to accretion rate changes. However, emission lines can be originated from other processes apart from accretion and thus their variability traces additional mechanisms and timescales. A good example is the bulk of the [OI]6300 emission line, which in HAeBes mainly probes to the low-density, surface layers of the protoplanetary disks [113] although its luminosity still correlates with L acc and serves to provide a rough estimate of the accretion rate [61]. In fact, the simultaneous monitoring of the near-UV excess-directly related to accretion-and several line luminosities show that they can evolve differently in HAeBes (e.g., [61,62]). The careful monitoring of specific CTTs actually reveals significant time delays between the moment when the material shocks onto the stellar surface and the time when such an accretion event is reflected by several emission lines at different wavelengths [114].
On the other hand, the previously discussed direct and indirect methods to estimate accretion rates in HAeBes can be considered an extension from similar methodologies previously applied to CTT stars. However, it is worth mentioning that not all methods used in the low-mass regime can be extended to the HAeBes, which is summarized as follows: • Spectroscopic line veiling: CTTs show optical photospheric absorption lines smaller in depth than observed in WTTs or low-mass stars in the MS, which can be explained from the contribution of the hot accretion shock. Indeed, by removing the contribution of the stellar photosphere to the observed photospheric lines one can directly infer a value ofṀ acc (see e.g., [115][116][117][118] and references therein). In contrast, optical spectroscopic line veiling is not commonly observed in HAeBes, not because they are not accreting but due to the fact that the temperature of the accretion shocks is comparable or smaller to that of the stellar photosphere (∼10000 K) and therefore the contrast effect is negligible (see [56] for more details). Thus, optical spectroscopic line veiling is not a method to routinely derive accretion rates of HAeBes excepts perhaps for the coldest sources (see an example in [119]).
• Hα line width: The Hα line width at 10% of peak emission, W 10 (H α ), not only serves as a qualitative indicator of accretion in CTTs and brown dwarfs [117] but also as a rough quantitative estimator through an empirical correlation withṀ acc [39]. Based on the study of a HBe source, Boley et al. (2009) [120] already suggested that such an empirical correlation may not extend to higher masses, which was later confirmed by   [61] from the study of a larger sample of HAeBe stars. This work showed that while their typically large vsin i values are reflected by Hα emission broadening in possible agreement with MA, the influence of rotation is that important that theṀ acc -W 10 (H α ) correlation breaks and thus cannot be used for the HAeBe regime.
In summary, most-but not all-methods that are commonly used to derive accretion rates in the low-mass regime and are based on MA can be extended to the HAeBes with similar accuracies.

Boundary Layer Measurements of HAeBe Stars
As we have seen almost all efforts devoted to derive accretion rates in HAeBe stars rely on the MA paradigm, although the first estimates were actually based on the BL scenario and the observed near-IR excesses (see Section 2 and [47]). However, such a wavelength range does not reflect accretion in most sources (see e.g., the discussions in [56,121] and references therein). Alternatively, Blondel & Tjin A Djie (2006) [48] derived accretion rates of a relatively wide sample of HAeBes based on the observed UV excesses and the BL scenario. In this work, low resolution spectra from the International Ultraviolet Explorer (IUE) combined with optical information of dozens of late-type HAeBes (HAes and IMTTs) were reproduced in terms of the photosphere of a central star, an optically thick accretion disk, and a hot BL in-between the previous. As a result of such a modeling, stellar parameters and accretion rates were derived for most of the stars in that sample, constituting a unique database that can be compared to MA estimates mentioned above.
A direct comparison between the accretion rates derived by Blondel & Tjin A Djie (2006) [48] from BL and the ones inferred from MA shock modeling by  and Fairlamb et al. (2015) [61,63] suggests that MA estimates tend to be smaller than from BL for many stars, in agreement with a similar comparison carried out for TTs [23]. However, the stellar parameters and distances used in the BL and MA works differ-sometimes significantly-which could affect the previous conclusion. An alternative approach to compare accretion rates from MA and BL, avoiding the influence of stellar parameters and distances, is discussed next.  Table 2 in [48]) versus the L * values derived in that same work. In the BL formalism, L acc ∼ (1/2) × GM * Ṁacc /R * (see Section 6). TheṀ acc values used in this formula are the ones tabulated in Blondel & Tjin A Djie (2006) [48] with inclinations to the line of sight different than 90 • , when available. The edge-on values from that paper are taken otherwise, which can be considered to be lower limits (triangles in Figure 3). In addition, for the stars with different estimates the averages are computed, the errorbars indicating the largest difference with respect to the individual values. The fit to the L acc -L * values of HAes when accretion is estimated from the BL view in Blondel & Tjin A Djie (2006) [48] is indicated with a black line and has a slope ∼2. In contrast, the L acc -L * trend followed by HAes when their accretion luminosities are inferred from MA is overplotted in green. This is equal within errorbars to the one followed by TTs and has a smaller slope ∼1 (see the bottom right panel of Figure 2 and [65,66] for details).  [48] were used to infer the L acc values, as well as their data for M * and R * to derive the final MA-based values foṙ M acc = L acc R * /GM * . The comparison between BL-and MA-based estimates ofṀ acc shows that the former tend to be larger for stars accreting at relatively high rates, which constitute ∼ 41% of the stars in the sample analyzed. On the opposite, the weak accretors showṀ acc (BL) <Ṁ acc (MA), representing ∼ 21% of the stars in the sample. The rest of the stars, ∼ 31%, show accretion rates that do not differ significantly regardless of the accretion scenario. The mismatch between MA and BL estimates is in last term due to differences between the underlying assumptions concerning energy transformation under both approaches (e.g., [66]). In MA the gravitational energy is carried through the accretion funnels and released at hot spots in the stellar surface, whereas in BL the kinetic energy of the rotating disk heats the BL when material is drastically decelerated in this region. Further details about BL in comparison to MA are discussed in the next section.

The Ultraviolet Link
The analysis in the previous section illustrates the importance of adopting a given accretion paradigm to interpret the data, as different scenarios lead to different values of the mass accretion rates. In turn, this has a potential effect on the shape of the L acc -L * correlation (or the roughly equivalent one betweeṅ M acc and M * ) and other related scaling relations ( Figure 3 left; [65,[122][123][124], and references therein), the value of accretion-based disk masses [69,121,[125][126][127], the disk dissipation timescales and driving processes including planet formation (e.g., [9,121,[128][129][130][131], and references therein), the removal of angular momentum (e.g., [132,133], and references therein), or even on the estimates of SFRs characterizing large star-forming regions when inferred from the sum of individual stellar accretion rates [7,8], to cite some. Although so far the bulk of the evidence supports MA as the driving mechanism in HAe stars, the possibility that the alternative BL scenario plays a role in these objects, and especially in HBe stars, needs to be further explored. In this section, simple MA and BL shock models are compared to argue that the study of the UV wavelength region could contribute to discriminate between both competing accretion scenarios.
The BL is geometrically described as an annulus around the star with radial thickness δ. Following Lynden-Bell & Pringle (1974) [14] δ ∝ ν 2/3 GM * The stellar and Keplerian angular velocities (at the star's surface, i.e., the break-up velocity; Ω * /Ω K < 1) are Ω * = vsin i/R * sin i and Ω K = (GM * /R 3 * ) 1/2 , respectively. The right term in Equation (2) is derived assuming that the disk viscosity can be expressed as ν ∝ (GM * R * ) 1/2 [14,134], and it reflects the fact that the radial thickness of the BL increases for higher stellar to Keplerian velocity ratios. The exact value of δ 0 (and therefore of δ) is difficult to constrain observationally. Following Blondel & Tjin A Djie (2006) [48], it will be assumed a typical (median) value δ 0 = 0.03 by default, although this is a free input parameter in principle ranging from 0.01 to 0.5 (see e.g., [17]).
Concerning the energy balance, it is assumed that the kinetic energy of the disk material rotating Keplerian is transformed to heat the BL [14,66], which in zeroth order approach is represented by an optically thick region characterized by a single blackbody temperature T BL and emitting area 4πR * δ (see [48,135]). Please note that all the energy is assumed to be released in the BL, and the part devoted to spin-up the star is neglected. This is observationally justified considering that the majority of A and B stars rotate much slower than the break-up velocities at least at the beginning of the MS [136,137].
The total (non-extinct) flux per wavelength unit coming from the star-disk system (F λ ) can be divided into three components: the flux from the boundary layer (F BL λ ), the photosphere of the star (F * λ ) and the disk (F disk λ ) Given that the contribution of the disk in the UV is negligible, it will not be considered here anymore. The first term in the previous equation can be expressed as (see [15,134]), where B λ (T BL ) is the blackbody radiance characterizing the boundary layer. The factor (π + 2γ)cos i describes the effect of partial screening of the boundary layer by the stellar surface due to inclination of the system in the line of sight, with γ derived by geometrical arguments as (see [134]). Therefore, the BL flux received increases for larger radial thicknesses and for inclinations closer to pole-on, when the BL emitting region is less occulted by the stellar surface. Regarding the second term in Equation (4), this is where the photospheric flux is corrected by the fraction of the stellar surface that is covered by the BL annulus. This is supposed to have a physical height similar to its radial thickness δ [15,48]. Therefore, the photospheric flux received decreases for inclinations closer to edge-on, when the stellar surface is more occulted by the BL. For a given star with stellar parameters M * , R * (or log g), projected rotational velocity vsin i and inclination to the line of sight i, the thickness of the boundary layer δ is estimated from Equation (2). From this, different mass accretion rates provide different boundary layer temperatures through Equation (3). Then the flux contribution of the boundary layer and that of the stellar surface (using a Kurucz template with a given stellar temperature and surface gravity) is derived through Equations (5) and (7), respectively. Finally, the expected total flux is obtained from Equation (4). As mentioned above, here we will focus on the UV range, for which the disk emission in the previous equation can be neglected.
Concerning MA shock modeling, the same prescriptions described in  [61] will be used here in order to compare the MA and BL predictions in the UV. There are two major differences with respect to the BL model described above. First, in the BL model the emitting region due to accretion is an equatorial annulus with thickness δ around the star, and the received flux depends on the inclination of the source. In contrast, in MA the emitting region is described from the accretion spots assumed to be located at high latitudes (i.e., inclination independent), and the corresponding fraction of the stellar surface covered by such spots ("filling factor"). Secondly, for a fixed value ofṀ acc and a given set of stellar parameters, the free parameters in the MA model are the inward flux of energy carried by the accretion columns, (F , normally expected to be around 10 12 erg cm −2 s −1 ; see [56]) and the disk truncation radius (R i , that depends on the specific star but it is typically ∼2.5 R * for HAeBes, see Section 4). These parameters determine the temperature and filling factor characterizing the accretion columns and shocks (larger values of F and R i provide larger values of T col and smaller filling factors, respectively). In turn, in the BL model the free parameter is the size of the emitting region inferred from δ 0 , which also determines the BL temperature (larger values of δ 0 provide smaller values of T BL and vice versa). The interested reader can consult   [61] for more details on the MA shock model.
Here, the MA and BL predictions in two regions of the UV spectra will be compared to each other. The excess in the Balmer region is defined as where F U /F phot U represents the ratio between the mean total, accretion-contributed flux and the mean photospheric flux in the wavelength region 3500-3700 Å (i.e., the U photometric band), where both are normalized at around 6000 Å (V-photometric band). Kurucz templates are again used to represent F phot U . Analogously, to quantify the excess at a shorter wavelength it is defined with F 2125 /F phot 2125 being the ratio between the mean total, accretion-contributed flux and the mean photospheric flux in the wavelength region 2000-2250 Å. Such a wavelength region has been chosen arbitrarily only for comparison purposes with ∆D B , but a different UV range could be analyzed instead. It should be noted that the effect of extinction is not considered either in the previous definitions or for the subsequent discussion based on analytical expressions. In addition, and assuming 10% and 5% errors in the total and photospheric fluxes, a typical uncertainty for the above defined excesses is ∼0.1 magnitudes, which represent the maximum errorbars that could be derived from broadband photometry. Such errors could be significantly reduced if e.g., high precision photometry or spectra are used. Figure 4 compares the excesses predicted by the MA and BL shock models applied to three "typical" HAeBes representing a late-type HAe star (top panels), and early-type HAe star (mid-panels), and a mid-type HBe star (bottom panels). Two different values for the gravity at the stellar surface are considered; log g = 4 and 3 for the left and right panels, respectively. Aiming to extract some general conclusions, typical values of F = 10 12 erg cm −2 s −1 and R i = 2.5R * (for MA), and δ 0 = 0.03, vsin i = 100 km s −1 , and i = 45 • (for BL) have been adopted in all cases. The analysis of Figure 4 provides the following general conclusions: • The relation between the excesses and the accretion rate is strongly dependent on the stellar properties not only for MA (see also Figures 1 and 9 in [61,63] respectively) but also for the BL models. A given excess can correspond to accretion rates different by orders of magnitude, depending on the stellar parameters (T * , M * /R * ) of the source. Additional dependencies in the BL model are discussed below.

•
The differences between predictions from MA and BL become more significant (i.e., above errorbars from broadband photometry) for high excesses/accretion rates. Those differences are generally larger for the Balmer excess than for the UV excess, at least for the space parameters explored here. For instance, the Balmer excess predicted by MA in the mid-left panel of Figure 4 is ∼3 times larger than from BL, forṀ acc ∼ 3 × 10 −6 M yr −1 .

•
For the general case analyzed here BL requires higher accretion rates than MA to reproduce a given, large enough Balmer excess. In other words, for relatively large Balmer excesses accretion rates predicted from MA are lower limits to the corresponding from BL, in agreement with the discussion in Section 5.

•
The ratio between the UV and Balmer excesses (∆ UV /∆D B ) predicted by both models tends to differ significantly in most cases. For instance, the ratio predicted in the mid-left panel of Figure 4 is close to unity from MA, while it can reach a factor ∼2 from BL. Similarly, for a given star ∆ UV /∆D B can be < 0 from MA, and > 0 from BL (bottom panels in Figure 4). This particular difference between predictions represents a good opportunity to observationally compare the two competing models.
Although the previous analysis indicates that the competing MA and BL scenarios could be eventually distinguished from simultaneous observations at different UV wavelengths, such an analysis does not consider the specific properties of each star-disk system, which, as mentioned before, are crucial to properly reproduce the observations of a given object. Concerning the BL modeling in particular, Figure 5 shows the influence of varying the size of the BL (left) and the inclination to the line of sight (right). While different inclinations provide excesses that are roughly consistent to each other at least within errors from broadband photometry, the value assumed for the BL size is more critical, and changing it by 1 dex results in significant variations above photometric errorbars in the expected excesses both in the Balmer region and at shorter wavelengths. It is recalled that the reason why the large Balmer excesses shown by several HBe stars cannot be explained from MA is that the corresponding filling factors representing the fraction of the stellar surface covered by accretion spots should be above 100% [61,63], which obviously is not physically possible. In contrast, the large Balmer excesses shown by these HBe stars could be potentially reproduced from BL by changing δ 0 (and thus T BL ) somewhat arbitrarily, given the lack of observational or theoretical constraints on the BL size. However, as has been shown above, the simultaneous measurements of continuum excesses at different wavelengths limits the possibilities that could be reproduced from BL, potentially constituting a more stringent test of this scenario.
On top of the previous caveats concerning the specifics of each star-disk system, future analysis in the UV must also carefully address the problem of extinction correction. This has not been considered in the previous analysis, but it is critical to properly deal with observational data in that wavelength region, as extinction fundamentally affects the excesses that want to be reproduced.
The Hubble Space Telescope's Ultraviolet Legacy Library of Young Stars as Essential Standards (ULLYSES) 5 will soon observe a large sample of young stars in the UV probing the lowest and the highest stellar temperatures (K-M and early O-B stars). Although ULLYSES will not cover the intermediate-mass regime, the resulting dataset in combination with related ground-based efforts will surely provide valuable strategies to best deal with related archival data of HAeBe stars. In the longer term, >100 h with the future World Space Observatory 6 will in principle be devoted to observe HAeBes in the UV and tackle the problem of the best suitable accretion scenario in these stars.

Concluding Remarks
In this work our estimates of the mass accretion rates in HAeBe stars have been reviewed. Guided by previous works devoted to the better known low-mass regime, we have seen that MA also seems to be a valid scenario that explains the properties of the late-type HAes and is capable of providing accretion rates ∼ an order of magnitude larger than for ∼1 M CTTs. Still, future work is necessary to reach a consensus on which accretion scenario is best suited for the HAeBe regime. Two main lines of research have been suggested to clarify this open issue. First, it is necessary to reduce the detection limits involved in the determination of magnetic fields in HAeBes, which currently are of the same order than the minimum magnetic fields required to drive accretion magnetospherically in these sources. Secondly, more efforts should be devoted to the theoretical and observational developments of alternative accretion scenarios. Although it has been shown from simple accretion shock models that the UV region could help to disentangle between MA and the competing BL scenario, proper modeling of the BL region is still pending. Indeed, radiative transfer models able to predict the shape and strength of emission lines from BL would be particularly useful. Alternatively, other possible accretion paradigms different than MA and BL could work in HAeBes (see e.g., [138]).
Even assuming that MA can be extended from the CTTs to the late-type HAes, MA shock modeling is not capable of reproducing the Balmer excesses of several early-type HBe stars, which require a different approach [61,63]. If accretion in these sources occurs through a BL, the analysis presented in this work suggests that they would be accreting at very high rates, larger than inferred from the MA-based correlations with the emission lines or the stellar luminosities. For instance, the large photometric Balmer excesses of HBe stars that cannot be explained from MA are >0.5 magnitudes and could reach >1 magnitudes for some sources [61], in principle corresponding to rough accretion rates >10 −3 M yr −1 using large enough BL sizes ( Figure 5, left). Such huge accretion rates are indeed higher than in principle expected for HBe sources (e.g., [10], and references therein), which by definition have already dissipated their natal envelopes and are very close to the MS. Although this review has been focused on accretion, alternative physical processes that do not invoke infalling material may also be capable of explaining the observed excesses in HBes, and perhaps several other observational differences between these and lower-mass stars. In particular, there is evidence that photoevaporation plays a more important role in HBes than in HAes (see e.g., [69,100,139], and references therein). In turn, photoevaporated disks should show relatively small accretion rates according to models (e.g., [140,141], and references therein), again suggesting that > 10 −3 M yr −1 may be unrealistically large for HBes. Therefore, it is worth exploring the possibility that the observed Balmer and UV excesses in HBes are not explained mainly in terms of accreting material shocking onto the star, as for HAes and CTTs, but perhaps in terms of photoevaporative outflows shocking onto the disks.