Observational constraints on the origin of the elements II. 3D non-LTE formation of Ba ii lines in the solar atmosphere

Context. The pursuit of more realistic spectroscopic modelling and consistent abundances has led us to begin a new series of papers designed to improve current solar and stellar abundances of various atomic species. To achieve this, we have began updating the three-dimensional (3D) non-local thermodynamic equilibrium (non-LTE) radiative transfer code, Multi3D, and the equivalent one-dimensional (1D) non-LTE radiative transfer code, MULTI. Aims. We examine our improvements to these codes by redetermining the solar barium abundance. Barium was chosen for this test as it is an important diagnostic element of the s-process in the context of galactic chemical evolution. New Ba II + H collisional data for excitation and charge exchange reactions computed from first principles had recently become available and were included in the model atom. The atom also includes the effects of isotopic line shifts and hyperfine splitting. Method. A grid of 1D LTE barium lines were constructed with MULTI and fit to the four Ba II lines available to us in the optical region of the solar spectrum. Abundance corrections were then determined in 1D non-LTE, 3D LTE, and 3D non-LTE. A new 3D non-LTE solar barium abundance was computed from these corrections. Results. We present for the first time the full 3D non-LTE barium abundance of $A({\rm Ba})=2.27\pm0.02\pm0.01$, which was derived from four individual fully consistent barium lines. Errors here represent the systematic and random errors, respectively.


Introduction
Barium is key element that is used in heavy element studies in stars. Its abundance patterns in the halo, in field stars, and in clusters have been carefully measured over the past several decades. Barium, like most other heavy elements, mostly forms via a series of neutron captures through either the rapid (r-) process or slow (s-) process channels. These two neutron capture channels have very different sites. After the discovery and analysis of 2017gfo (the electromagnetic counterpart of GW170817 Valenti et al. 2017), it is highly probable that the r-process mostly occurs in neutron star mergers (Thielemann et al. 2011). Conversely, the majority of barium in the Sun (81% Arlandini et al. 1999) ostensibly formed via the s-process in thermally-pulsing asymptotic giant branch (TP-AGB) stars (Smith & Lambert 1988). However, other sites for the s-process and r-process do and most likely exist. Naturally, the barium isotope ratio, f odd 1 , of a star is a useful quantity as it provides precise information on the s-and r-process contribution, but is exceedingly difficult to measure. Therefore, this information is only measured in some thick disk and halo stars, where this parameter is most interesting (Magain & Zhao 1 f odd ≡ N( 135 Ba) + N( 137 Ba) /N(Ba) 1993; Magain 1995;Mashonkina et al. 1999;Gallagher et al. 2010Gallagher et al. , 2012Gallagher et al. , 2015.
Most abundances, save those such as lithium that are measured in absolute units, are measured relative to the solar abundances. This helps to mitigate systematic errors within spectroscopic abundances and yields extra information about stellar populations, evolutionary stages and ages, that measurements in absolute units might not. As a result, the solar abundances are extremely important to stellar astrophysics. Consequently, very accurate measurements of the solar abundance are needed, which employ sophisticated model atmospheres and spectrum synthesis techniques such as 3D hydrodynamics and non-local thermodynamic equilibrium (non-LTE) physics. In recent years, with the development of faster and larger computers, it has been possible to develop and implement these methods (Asplund et al. 2003;Steffen et al. 2015;Klevas et al. 2016;Amarsi et al. 2016;Mott et al. 2017;Amarsi & Asplund 2017;Nordlander et al. 2017).
One of the main aims of this paper series is to report on our development of the one dimensional (1D) and three-dimensional (3D) statistical equilibrium codes -MULTI 2.3 (Carlsson 1986) and MULTI3D (Leenaarts & Carlsson 2009) -as we include new or better physics into their program flows. Given how impor-Article number, page 1 of 11 arXiv:1910.03898v3 [astro-ph.SR] 6 Feb 2020 A&A proofs: manuscript no. Ba_nlte tant barium is to galactic chemical evolution studies because it traces the impact of neutron-capture nucleosynthesis, we present a thorough analysis of the solar barium abundance using a handful of Ba ii optical lines computed using the two statistical equilibrium codes and the same barium model atom.
The statistical equilibrium of Ba ii has already been a subject of several detailed studies (Mashonkina & Bikmaev 1996;Mashonkina et al. 1999;Shchukina et al. 2009;Andrievsky et al. 2009;Korotin et al. 2015). The first such study was conducted by Gigas (1988) in Vega. There are, however, important differences between our work and these earlier studies. First, we use the new quantum-mechanical rates for transitions caused by inelastic collisions with hydrogen atoms from Belyaev & Yakovleva (2018). We also examine the impact dynamical gas flows have on Ba ii by utilising a 3D radiative hydrodynamical model to compute full 3D non-LTE radiative transfer, as well as 3D LTE, 1D LTE, and 1D non-LTE. Ab initio collisional damping from Barklem et al. (2000) was included in the linelist.
It has been observationally confirmed that the Ba ii resonance line at 4554 Å is sensitive to the chromospheric effects 2 , and so naturally a polarised spectrum of the resonance line is also sensitive to the quantum interferences (see, e.g. Kostik et al. 2009;Shchukina et al. 2009;Belluzzi & Trujillo Bueno 2013;Smitha et al. 2013;Kobanov et al. 2016), however, this is beyond the scope of this paper.
The paper is structured as follows. In Sect. 2 we describe the observations, we detail the model atmospheres, model atoms and spectral synthesis codes; in Sect. 3.2 we discuss the impact that various model assumptions have on our results; in Sect. 5 we describe the analysis and results from our Ba ii line analysis; and in Sect. 6 we summarise the study.

Solar spectrum
The solar spectrum is taken from the Kitt Peak National Observatory (KPNO) solar atlas published by Kurucz et al. (1984). This solar atlas covers the spectral range of 3 000 to 13 000 Å at a typical resolution R ≡ λ ∆λ = 400 000. Although newer solar spectra exist such as the PEPSI spectrum provided by Strassmeier et al. (2018), we chose to work with the former atlas as it has a very high resolution, roughly twice that of the latter. Nevertheless, comparisons of these two spectra have previously been made and they were found to be in very good agreement with one-another (Osorio et al. 2019).

1D model atmosphere
We use the MARCS model atmosphere that was computed for the Sun from the opacity sampled grid published in Gustafsson et al. (2008). The solar parameters of this model are T eff / log g/[Fe/H] = 5777/4.44/0.00 and include a mixing length parameter, α MLT = 1.50. The solar composition used to compute the model opacities are based on those published in Grevesse et al. (2007).

3D model atmosphere
For the work presented in this study we make use of the solar stagger (Nordlund et al. 1994;Nordlund & Galsgaard  1995) model with stellar parameters T eff / log g/[Fe/H] = 5777 K/4.44/0.0, from the stagger model atmosphere grid (Collet et al. 2011;Magic et al. 2013). A 3D model consists of a series of computational boxes that represent a time series, which are commonly referred to a snapshots. These snapshots are selected from a larger time series of snapshots that are produced from the stagger code and are selected at a time when the simulation has reached dynamical and thermal relaxation. For our purposes -and for the sake of time -we have chosen to work with five snapshots, each consisting of 240×240×230 grid points which cover a geometrical volume of 7.96 × 7.96 × 3.65 Mm in the x, y, and z dimensions, respectively. These are used as independent input for our statistical equilibrium code MULTI3D (Sect. 2.5) and then the output are averaged together, therefore applying the ergodic approximation that averaging in time is equivalent to averaging over space. In this case it is assumed that averaging in time is equivalent to averaging across the disc of the Sun. The emergent fluxes from these snapshots have an equivalent width variance of only ∼ 0.75 mÅ, suggesting that including further snapshots to the study will not greatly improve the results presented in this study, only increase the computational times.
Line opacities were collected from the MARCS database and are sorted into 12 opacity bins. Continuous absorption and scattering coefficients are taken from Hayek et al. (2010). Importantly, and unlike an equivalent 1D model, 3D models provide x, y, and z velocity fields for every voxel meaning that postprocessing spectrum synthesis codes provide more accurate approximations for the Doppler broadening, including asymmetric line profiles, which result from these gas flows.
We also make use of the averaged 3D model to help make qualitative comparisons between the full 3D and 1D models, however, we do not use it with MULTI 2.3 or MULTI3D. A 3D model is computed from a 3D model by spatially averaging the thermal structure of the 3D computational box over surfaces of equal Rosseland optical depth. As this can be performed in several different ways, comparing results from different studies that do not specify their averaging techniques is ultimately selfdefeating. Figure 1 depicts the 3D solar temperature structure (blue 2D histogram), along with the 1D MARCS (dashed line) and 3D (solid red line) temperature structures. It is clear that the average temperature of the full 3D model and the 1D model are fairly consistent in the outermost regions of the atmosphere (as seen by comparing the 1D with the 3D model). However, in deeper regions of the models -where the continuum usually forms (log τ ROSS ≈ 0) -the models begin to diverge. This is mostly due to the differences between the convection indicative in the 3D hydrodynamic model atmosphere -which the 3D model traces -and the treatment of convection theory (in this case the mixing length theory) in the 1D model atmosphere.

MULTI 2.3
MULTI solves the equations of radiative transfer and statistical equilibrium in 1D geometry with 1D model atmospheres. The latest release of MULTI is MULTI 2.3. However, we have made several minor changes to MULTI 2.3 for our purposes including, the ability to compute the detailed balance for charge transfer processes between ions and hydrogen. We include a fixed microturbulence value of 1 km s −1 in our computations with the solar MARCS model. The flux data were computed using five µ-angles, assuming a Gaussian quadrature scheme taken from Lowan et al. (1942).

MULTI3D
MULTI3D is an message passing interface (MPI)-parallelised, domain-decomposed 3D non-LTE radiative transfer code that solves the equations of radiative transfer using the Multi-level Accelerated Lambda Iteration (ALI) method (Rybicki & Hummer 1991 for 3D model atmospheres. Every element that is modelled by MULTI3D is assumed to have no effect on the model atmosphere, as it is in MULTI 2.3. This is a good assumption for barium as it is not an electron donor nor does it have a high impact on the overall opacity, unlike magnesium or iron, for example. At present, it will accept three types of 3D model atmospheres formats as direct input, including those computed using Bifrost (Gudiksen et al. 2011), and stagger. While Bifrost models are read using MPI IO, the stagger models are, at present, not read this way due to complications in converting byte ordering. However, the added delay to the code's run time is minimal, and only becomes noticeable when MULTI3D is run on several hundred CPUs. In addition to these two types of model atmosphere, the code will also accept any 3D model formatted so that the temperature, T , density, ρ, electron number density, n e , and x, y, and z velocity fields are supplied on a Cartesian grid that is both horizontally periodic and equidistantly spaced. Therefore, it is relatively straightforward to convert almost any 3D model to this input format for MULTI3D.
We have introduced new coding for computing fluxes inside MULTI3D along with the appropriate post-processing routines designed to extract the flux data. All of the output flux data computed for the work presented here was calculated using a Lobatto quadrature scheme and the appropriate corresponding weights (Abramowitz & Stegun 1972). At a later stage of this paper series, other quadrature schemes will be introduced, as well as internal routines that will compute fluxes inside MULTI3D and write them as output.
MULTI3D is now capable of accepting model atoms that include hyperfine structure (HFS) and isotope shift information for any atomic transition. This means that lines with highly asymmetric profiles, caused by these effects, can now be adequately modelled by MULTI3D. To test this upgrade, and to test that we could limit the impact of systematic errors dominating the abundances and abundance corrections we provide, we compared 1D spectra computed by both MULTI 2.3 and MULTI3D. This was conducted only for the vertical intensity (µ = 1), using a small test barium model atom, under the assumption of LTE. We use the same opacity sources, and same input model atmosphere. Systematic differences in the equivalent widths of less than 2.2 % were found between intensities computed with MULTI 2.3 and intensities computed with MULTI3D. This translates to abundances differences much less than 0.01 dex. The reason for these small differences is likely because of the way each code solves the radiative transfer equation; MULTI3D uses a direct 1D integration of the radiative transfer equation when computing spectra from 1D model atmospheres 3 , while MULTI 2.3 utilises a faster Feautrier method. However, abundance uncertainties found here are far smaller than the errors we report in Sect. 5. Therefore, we were satisfied that comparing 1D output from MULTI 2.3 with 3D output from MULTI3D was adequate.
We ran MULTI3D in short characteristic 3D solver mode and used the solar stagger model as input. The stagger model's xy grid points were scaled down by a factor of 64 from 240 × 240 to 30 × 30 grid points using a simple bilinear interpolation scheme. Significant tests conducted in the first paper in this paper series, , henceforth, Paper I), revealed no significant loss of information in the horizontal gas flows that affected the line profiles in any noteworthy way. The horizontal components were also assumed to be periodic so that rays with very low µ angles could be computed without encountering a horizontal boundaries. The vertical grid size remained consistent with the original model atmosphere at 230 grid points.

Model of Barium atom
The model atom of barium is constructed as follows. The energy levels for the Ba i and Ba ii levels are extracted from the NIST database. Of these, we include eight energy states of Ba i up to the energy of 2.86 eV, and all available levels of Ba ii up to 9.98 eV. Fine structure is retained for the three lowest terms of Ba ii: 6s 2 S (ground state), 5d 2 D (∼ 0.65 eV), and 6p 2 P • (∼ 2.6 eV) (Table 1). Transitions between these terms are typically used in the barium abundance analysis of cool stars. Other levels are merged into terms, and their energy levels are represented by the weighted sum of the individual components (weighted by the statistical weights of the levels). In total, the model is comprised of 110 states and is restricted by the ground state of Ba iii at 15.2 eV for a total of 111 levels. Note that no lines of Ba i are observed in the spectra of FGK stars. For the Sun n Ba i /n Ba ii ≈ 10 −4 − 10 −6 (see Sect. 3.1), hence, Ba i is a minority species. As such, a detailed treatment of the neutral stage is of no importance to the statistical equilibrium Ba ii, which is the majority species in these stars. Figure 2 depicts the energy levels in Ba ii and transitions among them in the form of a Grotrian diagram. We have colour-coded the four lines we use here as follows: gold represents the Ba ii 4554 Å line; green represents the 5853 Å line; blue represents the 6141 Å line; and red represents the 6496 Å line. As there are too many energy levels in the model atom to accurately depict without overlapping energy states, this figure should be used for qualitative assessments only.
The radiative bound-bound transitions for Ba ii were extracted from the Kurucz database, 26.03.2017. We also compared the data with the NIST database. For the combined terms, the lines were merged and the transition probabilities coadded as described in Bergemann et al. (2012a). The oscillator strengths of the four diagnostic lines used for the present analy-A&A proofs: manuscript no. Ba_nlte Notes. The wavelengths are given in air. The equivalent widths correspond to the measurements in the solar KPNO flux atlas and are given in mÅ. The Van der Waals broadening parameters are taken from Barklem et al. (2000). The log g f values reported here are derived from the very accurate transition probabilities taken from De Munshi et al. (2015) and Dutta et al. (2016). In total, the model contains 284 spectral lines in the wavelength range from 1330 to 202930 Å. The transitions with log g f < −10 are not included. Most of these lines are represented by nine frequency points, except for the four diagnostic lines that we use in the abundance analysis, i.e the lines at 4554, 5853, 6141, and 6496 Å. These lines were represented with a profile containing 301 frequency points. Several UV lines are rather strong. To test whether nine frequencies were enough to describe those strong lines, we ran a test using MULTI 2.3. The departure coefficients from the model atom were compared with an atom that contained 100 frequency points for two strong UV transitions at 2304.247 Å and 2341.429 Å. It was found that these transitions have no effect on the populations of the four barium lines of interest. This justifies the number of frequency points chosen for transitions that were not of interest to us for this study. Damping by elastic collisions with hydrogen atoms are computed using the α and σ parameters from Barklem et al. (2000) where available. When this information was missing, we used the Unsöld approximation, which was scaled by 1.5. The wavelengths are taken from Karlsson & Litzén (1999). We found that for three lines (NB: not the line at 4554 Å) there was a systematic offset to the solar spectrum. Once the solar spectrum was corrected for gravitational redshift they did match the observed line positions, although the 4554 Å resonance line has a slightly different shift. Karlsson & Litzén underline that this line might be slightly shifted in their measurements, due to the isotopic mix they used and self-absorption in this strong line. This shift is however expected to be at most 1 mÅ to the blue (Litzén, private comm.) not sufficient to explain the remaining offset we observe. We discovered that the excess shift in this line was due to convective effects. We shifted the 4554 Å line by 2 mÅ in 1D to the blue to match the observed position, whereas we did not need to shift the 3D profile. We also introduce HFS and isotopic shifts. They were computed using the solar abundance ratios of the five barium isotopes, see Eugster et al. (1969). The odd barium isotopes have non-zero nuclear spins that causes hyperfine splitting of the levels. The magnetic dipole and electric quadrupole constants for the five relevant energy levels were taken from Silverans et al. (1986), and Villemoes et al. (1993). The isotopic shifts are provided by van Hove (1982) for the 5853 and 6141 Å lines, by Villemoes et al. (1993) for the 6496 Å line, and by Wendt et al. (1984)  The radiative bound-free data are computed using the standard hydrogenic approximation (Kramer's formula). This is appropriate since the first ionisation potential of Ba ii is at 10 eV, and the energy levels, which may contribute to radiative overionisation at the solar flux maximum, have very low population numbers. Also on the basis of earlier studies with strontium (Bergemann et al. 2012b), which has a similar atomic structure, we do not expect that photo-ionisation is a significant non-LTE effect. In fact, earlier studies of barium in non-LTE showed that Ba ii a majority ion and is collision-dominated (see Mashonkina et al. 1999;Gehren et al. 2001;Bergemann & Nordlander 2014). In such ions, the statistical equilibrium is established by a competition of collisional thermalisation, photon losses in strong lines, and over-recombination. This will be discussed in detail in the next section.
One of the new features of our atom, compared to earlier studies mentioned above, is the treatment of collisions. In particular, we include the new quantum-mechanical rate coefficients for the inelastic collisions between the Ba ii ions and hydrogen atoms by Belyaev & Yakovleva (2017. To the best of our knowledge, the first study of barium that employs these detailed quantum-mechanical data for collisions with hydrogen was recently published by Mashonkina & Belyaev (2019) for the purposes of treating isotopes. The present paper is the first time these hydrogen-collision data have been employed for full non-LTE modelling. The data are available for 686 processes, and represent collisional excitation and charge transfer reactions, i.e. Ba + H ↔ Ba + + H − . The rate coefficients are typically large and may exceed 10 −8 cm 3 s −1 in the temperature regime relevant to modelling the solar atmosphere. Note that here the process ionisation refers to the ion we are interested in and a free electron, i.e. ionisation: Ba ii + H → Ba iii + H + e, but the process ion-pair formation reads: Ba ii + H → Ba iii + H − , that is, Ba ii loses its outer electron and it is bound with H. Its inverse process, mutual neutralisation, is the process when Ba iii gains an electron from H − . The same is valid when Ba ii is replaced by Ba i and Ba iii is replaced by Ba ii. Charge transfer reactions do not lead to a free electron, which is the case that is usually modelled by a Drawin's formula (Drawin 1968(Drawin , 1969. Excitation and ionisation by collisions with free electrons are computed using the van Regemorter (1962) and Seaton (1962) formulae. A study of the impact of different collisional rates was presented in Andrievsky et al. (2009, Sect. 3.2). No differences between these classical recipes were found. An earlier version of this model atom was used in . We have since updated the oscillator strengths.

non-LTE effects
In terms of departures from LTE, Ba ii is a collision-dominated ion (see Gehren et al. 2001, and references therein). The ionisation potential is too high for photo-ionisation to play a significant role in the statistical equilibrium (SE) in FGK type stars. On the other hand, the term structure of the ion, with several very strong radiative transitions, favours strong effects caused by line scattering. In particular, there is radiative pumping at the frequencies of optically-thin line wings, τ wing < 1, as long as at the line centre τ core > 1. This mechanism acts predominantly at larger depths and leads to over-population of the upper levels of the transitions, at the expense of the lower states. On the other hand, strong downward cascades occur higher up in the atmosphere, where the strong line cores become optically thin, τ core < 1. This mechanism de-populates the upper levels via spontaneous de-excitations and this downward electron cascade causes over-population of the lower-lying energy states. As the statistical equilibrium of Ba ii has been extensively studied in the literature (Mashonkina & Bikmaev 1996;Mashonkina et al. 1999;Short & Hauschildt 2006;Mashonkina et al. 2008), hence in what follows we will only describe the main features of the non-LTE line formation and discuss the differences with the earlier studies.
Quantitatively, this behaviour can be visualised as plots of the departure coefficients b i 4 as a function of the continuum optical depth at 5000 Å, log τ 5000 . Figure 3 depicts the b i behaviour for the solar MARCS model atmosphere. To facilitate the comparison with the detailed study by Mashonkina et al. (1999), we have chosen the same axis range as that paper. Thick coloured curves correspond to 5 energy states, which are involved in the radiative transitions listed in Table 1: the ground state of Ba ii, 6s 2 S, and the low-excitation terms 5d 2 D and 6p 2 P • . Thin grey dotted curves show all other energy levels of Ba ii in the model atom. It is interesting that despite major differences in the model and numerical methods, including the properties of the atomic 4 The departure coefficient is defined as the ratio of atomic number density for a given energy level i computed in non-LTE to that of LTE, b i = n i,non−LTE n i,LTE . model, model atmosphere, and the SE code, the agreement between our results and that of Mashonkina et al. (1999, see their Fig. 2) is very good. In particular, the Ba ii ground state is entirely thermalised throughout the full optical depth range and develops a very modest over-population only at log τ 5000 < −4. The first excited metastable state 5d 2 D is also close to LTE due to strong collisional coupling with the Ba ii ground state, although minor departures in the atomic number densities, of the order of a few percent, are seen at log τ 5000 < −3 and higher up. The term 6p 2 P • at ∼ 2.6 eV shows stronger deviations from LTE already in the deep layers, log τ 5000 ∼ −2, where the non-LTE population of the level is only about 80% of the LTE value (that is equivalent to b i = 0.8). The pronounced depletion of the term is caused by photon losses in the lines, connecting the ground state and the lowest metastable state with 6p 2 P • . This under-population increases outwards as the lines become optically thin. All energy levels above 6p 2 P • are over-populated at log τ 5000 < 0, the effect that Mashonkina et al. (1999) attributed to radiative pumping.
It is interesting to briefly discuss the importance of microphysical processes in the SE of Ba ii. Andrievsky et al. (2009) suggest that photo-ionisation cross-sections is the main source of uncertainty in the SE of Ba ii. This is only true for Ba i, however, no lines of the neutral atom are observed in the optical or infra-red spectra of FGKM stars (Tandberg-Hanssen 1964). Over-ionisation of Ba i has no effect on the population of Ba ii, as the ratios of number densities of two ionisation stages are n BaI /n BaII ∼ 10 −4 at log τ 5000 ≈ 0, and drops to ∼ 10 −6 in the outermost atmospheric layers in the solar model. This ratio is even more extreme in the atmospheres of metal-poor stars. For example, for a model atmosphere of an RGB star with T eff = 4600 K, log g = 1.6, and [Fe/H] = −2.5, n BaI /n BaII ∼ 10 −5 at log τ 5000 ≈ 0, but approaches ∼ 10 −10 close to the outer boundary at log τ 5000 ≈ −5. On the other hand, photoionisation in Ba ii is not important. The ion has a very high ionisation potential, and its well-populated energy levels with low excitation potentials have ionisation thresholds in the far-UV, at λ < 1000 Å, where radiative flux in FGK-type stars is negligibly small. We also recomputed the departure coefficients assuming σ photo /100 and σ photo × 100. This is a very conservative estimate of uncertainty in the cross-sections, when comparing to a very similar atom, Sr, for which detailed quantum-mechanical cross-sections are available from (Bergemann et al. 2012a). It was found that only the non-LTE populations of the Ba iii ground state change, but none of the important Ba ii levels show any difference with respect to our reference model.
A more important ingredient for the SE of Ba ii seems to be the accuracy of the data for inelastic collisional processes, in particular, those between Ba ii and H i atoms. Short & Hauschildt (2006) suggest that the non-LTE line profiles are invariant to a factors of 0.1 − 10 changes in the rates of transitions caused by collisions (NB they used approximate analytical formulae to represent these data). This may hold for a limited range of stellar parameters. For example, in the case of the Sun, using the Drawin's recipe or QM data does not give significantly different results. However, it is known that metal-poor stars are sensitive to non-LTE effects (Bergemann & Nordlander 2014). Therefore, it would be reasonable to assume that they would also be sensitive to different collisional recipes. This is particularly true in hydrodynamic model atmospheres as the decoupled non-local radiation field leads to a cooling in the outer regions of the temperature structure, relative to the equivalent 1D model (see e.g. Gallagher et al. 2016, their Fig. 1). We intend to explore this in the near future.

Line formation in hydrostatic and inhomogeneous models
We begin with the analysis of LTE and non-LTE formation of Ba ii lines in the 1D hydrostatic solar model. We will then extend the analysis to radiative transfer with 3D inhomogeneous models.
As discussed in the previous section, the non-LTE effects in Ba ii are primarily dominated by line scattering. Consequently, deviations from LTE in the line source function are expected to be significant. Since the ratio of the departure coefficients for the lower i and upper j level of the transitions is below unity for the diagnostic Ba ii lines, b j /b i < 1, the ratio of source function to the Planck function (see Bergemann & Nordlander 2014, for the derivation), also drops below unity. In other words, the source function in the line is sub-Planckian and the non-LTE lines profiles shall come out stronger than the LTE lines. In some cases (e.g. Bergemann et al. 2012a, for Sr), this effect is modulated by the change of the line opacity. However, since κ ν ∼ b i , and the population numbers for the lower levels of all Ba ii lines are essentially thermal throughout the line formation depths, the line opacity is very close to its LTE value. This simple analytical picture is confirmed by comparing the LTE and non-LTE line profiles (Fig. 4). The non-LTE effects are small and amount to the abundance difference of −0.05 (4554 Å) to −0.1 dex (5853Å). The other two Ba ii lines show similar behaviour. Figure 4 demonstrates that the 3D profiles are asymmetric and also shifted blueward relative to the 1D profiles, which is expected. This is a natural result of the convective motions inside the 3D model, that the 1D model cannot replicate (Löhner-Böttcher et al. 2018;Stief et al. 2019). This is particularly obvious in the three subordinate lines, where the HFS has far less impact to the line shape than it does in the resonance line, where asymmetries are seen in both 1D and 3D. The 3D profiles, for a given barium abundance, are consistently weaker than their 1D counterparts, both in LTE and non-LTE. Therefore, positive abundance corrections are going to be needed to reproduce the same equivalent width, foreshadowing larger 3D LTE and non-LTE abundances over the 1D LTE counterpart. Like was shown in the hydrostatic case above, deviations from LTE should be significant because of line scattering.
The lines of Ba ii in the solar spectrum are strong, with EW's from 207 mÅ (4554 line) to 68 mÅ (5853 Å) line. As such, they are not only extremely affected by damping, but are also blended. The 4554 Å line is blended by a Fe ii feature close the line core, but has little impact on the line. The subordinate lines also show some blending. Korotin et al. (2015) avoided, in particular, the 4554 Å line in their analysis of FGK metal-poor stars in stars where [Fe/H] > −1.0. On the other hand, Grevesse et al. (2015) included the 4554, 5853, and 6496 Å line in their analysis of the Sun. They also used the non-LTE corrections for Ba ii lines. We explore the effects of line blending and damping in Sect. 5.
The current photospheric solar abundance of barium derived by Grevesse et al. (2015) -who applied non-LTE corrections to their 3D LTE abundance -is 2.25 ± 0.03 ± 0.07 (where errors represent the statistical and systematic errors, respectively). The abundance of barium in CI meteorites (Lodders 2003) is 2.19 ± 0.03. The main goal of this paper is to explore whether ab initio atomic data from physical experiments and detailed quantum-mechanical calculations are successfully able to describe the spectrum of the Sun. Accordingly, we now move on to report our 1D non-LTE, 3D LTE and 3D non-LTE corrections relative to the average barium 1D LTE abundance.

Computing abundance corrections
A grid of abundances for all four barium lines were computed using MULTI 2.3. These were fit to the solar spectrum using a χ 2 code that treats abundances, macroturbulent broadening and wavelength shifts as free parameters (see Sect. 2.6 for details of line shifts). The macroturbulences found ranged from 1.5 − 1.9 km s −1 . The code also normalises the fit to a local continuum for each line using two patches of spectrum either side of the line. We fixed the rotational broadening to v sin i = 1.6 km s −1 (Pavlenko et al. 2012). The 1D LTE we attain represent the best statistical fit from this χ 2 code.
We computed the four barium lines in 1D non-LTE, 3D LTE and 3D non-LTE using three abundances; A(Ba) = 2.17, 2.27 and 2.35 -covering typical abundances reported in the literature for the solar barium abundance. The abundance corrections tabulated in Table 2 were determined by fitting the grid of 1D LTE profiles to the 1D non-LTE, 3D LTE and non-LTE lines so that their equivalent widths matched. Unsurprisingly, it was found that the corrections from all three abundances were identical, hence corrections we provide are robust against the typical abundance range found by most studies on the solar barium abundance. No extra broadening was added to any of the depicted lines. The 3D profiles are shown to be shifted relative to the 1D profiles. This is a natural consequence of the convective shifts indicative to a dynamic model atmosphere.

Results
Our best fit profiles are compared with the solar flux profiles in Fig. 5. The best fit 3D lines are computed by MULTI3D using the corrections given in Table 2. The 3D non-LTE abundances are remarkably consistent, apart from the 6141 Å line, which is approximately 0.12 dex larger than the other three. The reasons for this will become apparent by the end of this section. Unlike in the 1D case, the best fit 3D LTE profiles show large deviations in the line cores, relative to their non-LTE counterparts, yet their equivalent widths remain very similar.
The lines and abundances discussed this far have assumed that the barium lines are unblended in the solar spectrum. In reality this is not the case. In fact, the abundances all barium lines (particularly the 6141 Å line) are dependent on line blending, as we now discuss.

Line blending corrections
The four barium lines we model here suffer from the effects of blending with other atomic and molecular species. Naturally, this impacts the abundances we derive when we assume that the line is clean. This is what was done when synthesising the lines with MULTI3D and MULTI 2.3, as these codes do not currently possess the capability to synthesise blends. To examine this, we used the VALD3 5 database together with the barium line information extracted from the model atom to create new line lists for the 5 http://vald.astro.uu.se/ 1D LTE spectrum synthesis code, MOOG 6 . This code was chosen as the interactive plotting tool makes recomputing synthetic spectra and fitting it with observed data very simple. Lines were computed with and without and the abundance of the clean barium line was adjusted until its line strength matched the blended line.
The Ba ii resonance line at 4554 Å is the strongest line presented here. Naturally, it would dominate most of the line depression at this spectral region. It is found that when blends are included, the 1D LTE barium abundance must be reduced by 0.01 dex. The 5853 Å line abundance was reduced by 0.03 dex. The 6141 Å line was found to be severally affected by blending, as the abundance had to be reduced by 0.16 dex. The 6496 Å barium abundance had to be reduced by 0.04 dex.
We will use these abundance corrections to determine the barium abundance in all four paradigms in Sect. 5.4. First, however, it is important to determine how the barium lines are affected by systematic uncertainties, as we now present.

Line damping uncertainties
We have previously mentioned that differences in the radiative transfer solvers used by MULTI 2.3 and MULTI3D lead to extremely small systematic uncertainties in the barium abundance. These are small (< 0.01 dex) enough to be dwarfed by uncertainties associated with the van der Waals broadening parameters, now discussed. The barium lines synthesised here have varying  Table 2. Abundances provided in the panels do not represent the final barium abundance, as line blends are excluded here. See Sect. 5 and Table 3 for details. line strengths. This means that they will react differently to uncertainties associated with the van der Waals line damping parameters, σ and α. If we conservatively assume that there is a 10% uncertainty associated with the damping parameters tabulated in Barklem et al. (2000) then we can see how this affects the abundance measured in each line and derive a systematic uncertainty for each line. This does not affect the abundance corrections, as these systematics will affect all syntheses equally and so they cancel out, but it is important to determine any uncertainty associated with the abundances determined from these corrections.
We examined how both the cross sections (σ) and velocity exponent (α) affect the line strength by varying both separately and computing new lines using MULTI 2.3. It was found that varying α by ±10% had no impact on the line strengths of any of the lines, and hence we cannot attribute any abundance uncertainty to this parameter. However, a 10% uncertainty in σ was found to affect the line strength in all four lines. The 4554 Å line is the strongest line measured (207 mÅ). It is therefore reasonable to assume that this line will be the most affected by this damping parameter uncertainty. It was found that the barium abundance determined from this line varies by ±0.03 dex. The 5853 Å line is the weakest line in this study (68 mÅ). As such the uncertainty we assign to the van der Waals parameter of ±10% leads to a change in abundance of only 0.01 dex. The 6141 Å line is the second strongest line measured (126 mÅ). This leads to abundance variations of ±0.03 dex. Finally, the 6496 Å line is also fairly strong in the solar spectrum (102 mÅ). The uncertainty we assign the damping parameter varies the barium abundance found in this line by ±0.02 dex. The list of associated abundance uncertainties can be seen in column four of Table 3.
Uncertainties in line blends are also of concern when computing abundances. No uncertainty information is given by VALD, so we again conservatively assume that the log g f values of these lines have a 10% uncertainty. The abundances attained from the 4554 Å line with and without line blending were virtually identical. As previously mentioned, this is because the resonance line dominates line depression around this spectral region. Accordingly, uncertainties in blended lines of ±10% do not affect the barium abundance. While the 5853 Å line is the weakest line analysed, it suffers the least from blending. As such, there is also no sensitivity in barium abundance found from varying the blended lines. The blends around the 6141 Å line have a large impact on the barium abundance. Naturally, uncertainties in log g f lead to an uncertainty of ±0.02 dex. Therefore, the inclusion of blend uncertainties increases systematic uncertainty of the 6141 Å from 0.03 dex to 0.04 dex. Finally, the blending uncertainties around the 6496 Å line were not found to influence the barium abundance. A break down of the associated abundance uncertainties can be seen in column five of Table 3.

Oscillator strength uncertainties
The oscillator strengths ( f -values) of the four diagnostic lines used in the present study are taken from De Munshi et al. (2015) and Dutta et al. (2016). They also provide unique errors associated to each transition probability, which we convert in to oscillator strength uncertainties. The transition probabilities are extremely accurate, so the resulting uncertainties are very small. When these are included in our calculations we find that the propagated abundance error associated with the 4554 Å line is ±0.00 dex. The weakest line in our sample (5853 Å) has a propagated abundance uncertainty of ±0.01 dex. The most blended line (6141 Å) is found to have a propagated abundance uncertainty of ±0.01 dex. Finally, the 6496 Å line has an associated abundance error of ±0.01. The associated abundance uncertainties for each line can be found in column six of Table 3.

Solar barium abundance
We now present the barium abundance in four paradigms: 3D non-LTE, 3D LTE, 1D non-LTE and 1D LTE. We correct the 1D LTE abundances using the corrections in Table 2 and Sect. 5.1, and weight them by their uncertainties listed in Table 3 using an inverse-variance weighted mean ( ω i X i , where ω i = 1 σ 2 ). We report for the first time a full 3D non-LTE solar barium abundance of 2.27 ± 0.02 ± 0.01, where errors given here are the systematic uncertainties and the random error determined as the standard deviation found in the line-to-line scatter of the 3D non-LTE abundances, which can also be seen in Table 3. This value is 0.08 dex larger than the meteoritic barium abundance of 2.19 ± 0.03 determined by Lodders (2003). We therefore find a photospheric abundance that is slightly larger than that reported from meteorites.
Using the same method just described, we find that A(Ba) = 2.31 ± 0.02 ± 0.03 in 3D LTE. Errors again represent the systematic and random errors, like above. When we use 1D model atmospheres and apply non-LTE physics to the post-process spectral synthesis we find that the 1D non-LTE barium abundance is A(Ba) = 2.11 ± 0.02 ± 0.05. Finally, when we derive the barium abundance using the classical 1D LTE approach we find that A(Ba) = 2.24 ± 0.02 ± 0.02. This abundance is similar to the 3D non-LTE abundance.

Conclusions
We have computed new values of the solar barium abundance based upon results computed using a new barium model atom that includes quantum mechanical inelastic collisional rate coefficients between Ba ii and hydrogen. We computed the barium abundance using a static 1D model atmosphere and provide 1D non-LTE, 3D LTE and 3D non-LTE corrections to this value in Table 2. Using these corrections we also present new solar photospheric barium abundances for the first time in full 3D non-LTE, as well as in 3D LTE, 1D non-LTE and 1D LTE.
The summary of this work is as follows (NB that all abundances given below are done so by adding the corrections to each 1D LTE abundance and then calculating the inverse-variance weighted mean as described in Sect. 5.4): -The 3D non-LTE barium abundance was found to be A(Ba) = 2.27 ± 0.02 ± 0.01, which is 4σ larger the meteoritic abundance published by Lodders (2003). This may suggest uncertainties in the atomic data and/or that further physics is still missing from our analyses, such as the treatment of magnetic fields. On the other hand, Ba isotopic abundance anomalies are well-documented in CI meteorites (e.g. McCulloch & Wasserburg 1978;Arnould et al. 2007) and it is not clear whether the meteoritic value can be indeed directly compared to the solar photospheric estimate. Nevertheless, this value represents the best photospheric solar barium abundance available from the current state-of-the-art in spectral modelling. As a result, it provides a remarkably consistent abundance for the four diagnostic lines. -The 3D LTE abundance was determined as A(Ba) = 2.31 ± 0.02 ± 0.03. This value is larger than the meteoritic value given in Lodders and the 3D non-LTE abundance we determine, and the individual abundances are not as consistent. -The 1D non-LTE abundance of A(Ba) = 2.11 ± 0.02 ± 0.05 suggests that the barium abundance is depleted by 0.16 dex in the solar atmosphere relative to the full 3D non-LTE. The abundance is also smaller than the meteoritic result reported in Lodders and the inconstencies between lines are larger than they are in 3D non-LTE.
The 3D non-LTE and 1D LTE abundances are very similar for barium in the Sun, but larger than that given in Lodders (2003). Conversely, the 3D LTE abundance suggests that barium abundance is even larger in the Sun, while the 1D non-LTE abundance suggests barium is slightly lower than the meteoritic value. It is clear then that the inclusion or removal of realistic treatments of line formation physics or convection has opposing effects on the barium abundance; the former strengthens the barium lines and so depletes the barium abundance, while the latter weakens the barium lines and hence a larger barium abundance is required. Therefore, the exclusion of both physical processes in the 1D LTE paradigm masks each effect, providing similar values in each line as the actual values found in 3D non-LTE. Conversely, in our work on manganese we found that the 3D and Non-LTE effects do not cancel out, but rather the effects of Non-LTE are amplified in 3D calculations with hydrodynamical models.
The previous set of transition probabilities reported by Davidson et al. (1992) led to abundances values that were, in general, less consistent than those reported here and had larger A&A proofs: manuscript no. Ba_nlte Notes. Column 2 tabulates the 1D LTE abundances found when the lines were assumed to be clean. They are consistent with the profiles depicted in Fig. 5. Column 3 presents the abundance correction associated with line blending. The final four columns present actual abundances of each line when the appropriate corrections were added to Cols. 2 and 3 from Table 2. The weighted averages calculated in Sect. 5.4 are computed using these abundances and the systematic errors tabulated in Col. 7 from the errors in BPO theory (Barklem et al. 2000), line blend uncertainties, and oscillator strengths, f , which are presented in Cols. 4, 5, and 6, respectively. uncertainties associated to them, with the 5853 Å line being most uncertain and most inconsistent with the other three diagnostic lines. The latest published transition probabilities in De Munshi et al. (2015) and Dutta et al. (2016) represent the most accurate transition parameters published. As such, the barium abundances we find from each diagnostic line used here are all in very good agreement (once the blending corrections are included).
We have presented this work as part of a larger series of papers designed to report on the development of the MULTI3D spectrum synthesis code. Up until now we have added new coding that allows it to read standard stagger model atmospheres; include the effect of charge transfer between hydrogen and ions; compute flux data based on hard-coded quadrature schemes; and compute multi-component transitions caused by HFS or isotope splittings. Further physics and mathematical techniques will be added as the project progresses that will be presented in future papers in this paper series. We also plan to extend our work on barium within this paper series to include several metal-poor benchmark stars, where the present work will be important to the relative abundances we report.
This following tables tabulate the complete hyperfine structure information of the four barium diagnostic lines used in this study.