7.6.4 Solver

The solver is the part of the code that fits synthetic LCs (generated by the EB simulator described in Section 7.6.3) to the observed data in order to determine the physical parameters of the system. This is accomplished in two steps. First, the observed LC is compared against a library of precomputed synthetic LCs to find one or more synthetic LCs that resemble the observed one. Then, the physical parameters associated with the best-matched synthetic LCs are used as starting values in a local optimisation procedure that converges to the parameters that best fit the observed LC. These two steps are described in more detail in what follows.

Global optimisation

The determination of EB physical parameters from a LC is a nonlinear least-squares problem. As such, the solver needs initial values for the parameters in order to search for the minimum. Furthermore, there are in general multiple minima of the $\chi^{2}$ function, and a local optimisation method (such as the Levenberg-Marquardt algorithm used in Section 7.6.4) cannot converge to the global minimum unless the initial parameter vector is located inside the convex region of the global minimum.

Traditionally, when fitting an EB LC “manually”, the task of finding a good initial parameter vector is the responsibility of the scientist, who will make a judgement based on the geometric characteristics of the LC as well as on any additional information that may be available on the EB in question. In the case of Gaia, the need for an automated processing of a large number of objects precludes manual intervention, thus mandating the use of an algorithm capable of finding a good initial parameter vector in an autonomous fashion. The selected approach involves comparing the observed LC against a reference library of synthetic LCs, uniformly sampled in phase at 100 points, in order to identify the top $N_{\mathrm{LC}}$ synthetic LCs that “best match” the observed LC.

The matching procedure involves phase-shifting, interpolating and scaling the synthetic LC until the $\chi^{2}$ difference between observed and synthetic LC is minimised. The parameter $N_{\mathrm{LC}}$ is configurable at runtime. A number of $N_{\mathrm{LC}}>1$ candidate LCs is desirable because (i) the parameter vectors of the LCs in the reference library are located at variable distances from the parameter vectors corresponding to their respective minima, and (ii) there can be multiple parameter vectors that produce synthetic light curves passing about equally well through the observed data points, either due to the sparse sampling of the observed LC or due to inherent degeneracies of the EB problem. As a result, the single best-matching synthetic LC may not lead to the best solution at the local optimisation stage. In practice, it was found that the improvement in local optimisation starts diminishing after $N_{\mathrm{LC}}\sim 5$; a value of 10 was used for the Gaia DR3 run.

The reader may wonder at this point whether some sort of parametric (e.g., polynomial) representation of the observed and reference LCs might have offered some advantages over a pointwise comparison of the LCs. Indeed, a parametric description could in principle extract the essential features of a LC while filtering out some of the noise in the data. While such an approach is not ruled out in a future DR, it would require a good sampling of critical LC segments that constrain the system parameters, such as during eclipse ingress and egress or the eclipse midpoint. This is something that the Gaia scanning law cannot guarantee, although phase coverage will of course improve as the data acquisition interval expands with time.

The reference library consists of $\sim\!\!1.6\times 10^{6}$ synthetic LCs precomputed with the EB simulator for a range of physical parameters. It contains $\sim\!\!60\%$ detached systems, $\sim\!\!20\%$ semi-detached systems, and $\sim\!\!20\%$ overcontact systems. All overcontact, almost all semi-detached systems, and approximately 1/3 of the detached systems are circular. The non-circular detached systems make $\sim\!\!40\%$ of the total number of systems, with $\sim\!\!24\%$ of the eccentricities in the [0.01, 0.3] range, $\sim\!\!10\%$ in the [0.3, 0.6] range, and $\sim\!\!6\%$ in the [0.6, 0.8] range. The sampling was uniform within each range. The argument of periastron, $\omega$, was also uniformly sampled within [$0^{\circ},360^{\circ}$).

The rotational parameters are set to $f_{1}=f_{2}=1$ in the case of circular systems, and computed by Equation 7.44 in the case of eccentric systems.

The inclination ($i$) was sampled in the range of $[75^{\circ},90^{\circ}]$ for detached systems, and $[40^{\circ},90^{\circ}]$ for semi-detached and overcontact systems.

The mass ratios ($q$), effective temperatures ($T_{\mathrm{eff}}$) and fill factors ($s$) were sampled in two distinct ways:

1. 1.

For the overcontact systems as well as for 1/3 of the detached and semi-detached systems, they were sampled randomly within reasonable ranges of these parameters. Specifically:

• $q\in[0.2,50]$ in all cases.

• $T_{\mathrm{eff}}\in[3\,000\,\mathrm{K},30\,000\,\mathrm{K}]$ determined separately for each component of the detached and semi-detached systems as well as for half the overcontact systems; for the remaining half of the overcontact systems, the temperature was set to be the same for both components, which is the usual common-envelope situation.

• $s\in[0.05,1]$ for both components of detached systems; $s_{1}\in[0.05,1]$ and $s_{2}\in[0.999999,1]$ for the components of semi-detached systems; and $s_{1}\in[1,1.4]$ for one component of overcontact systems, with the fill factor of the other component being determined via the common envelope constraint (cf. Section 7.6.3).

2. 2.

For 2/3 of the detached and semi-detached systems, these parameters were randomly chosen to correspond to the masses, radii and temperatures of main-sequence components.

Comparing an observed LC against all $\sim\!\!1.6\times 10^{6}$ synthetic LCs of the reference library is a relatively costly operation that scales linearly with the number of LCs in the library. Performance was improved, without significant loss in the quality of the solutions, through the use of a self-organising feature map (Kohonen 1982), an unsupervised machine-learning technique that implements an artificial neural network to perform a projection of the high-dimensional parameter space onto a two-dimensional grid. This process leads to the creation of a number of light curve templates, whereby each template represents a subset of the reference library LCs. In this way, the reference library can be split into a number of clusters containing “similar-looking” LCs. The self-organising map then groups similar clusters into neighbouring units such that increasing distance between clusters implies increasing “dissimilarity” between the LCs that comprise each cluster. This procedure constitutes the training of the neural network and is performed once ahead of the operational run.

In the global optimisation stage, the process of finding the best-matching LC is then split into two steps. First, the observed LC is compared against all template LCs to find the template that best resembles it. Then, the observed LC is compared against all the synthetic LCs that belong to the cluster of that particular template, as well as to neighbouring clusters within a radius defined by a runtime parameter. For the network used in the Gaia DR3 run, $9\,025$ clusters were defined, each containing an average of $\sim\!\!180$ LCs and no more than $\sim\!\!600$; the search radius was set to 2. Therefore, instead of comparing against $\sim\!\!1.6\times 10^{6}$ LCs, the observed LC was first compared against $9\,025$ template LCs, followed by 25 comparisons against an average of $\sim\!\!180$ LCs, for a total of $\sim\!\!13\,500$ LC comparisons on average, or a speedup by a factor of $\sim\!\!120$.

At the end of global optimisation, the following LC parameters are output for each of the $N_{\mathrm{LC}}$ light curves:

• time $t_{0}$ of eclipse of star 1 by star 2,

• mass ratio $q$ $(M_{2}/M_{1})$,

• fill factors $s_{1}$ and $s_{2}$,

• asynchronous rotation ratios $f_{1}$ and $f_{2}$ (different from 1 for eccentric systems only),

• inclination $i$,

• eccentricity $e$,

• argument of periastron $\omega$, and

• effective temperatures $T_{\mathrm{eff},1}$ and $T_{\mathrm{eff},2}$.

It should be noted that at this stage of the processing, the indices 1 and 2 of the stellar components are defined so that, at $t_{0}$, star 1 (the primary) is eclipsed by star 2 (the secondary). This is different from the definition in the final output, which is that component 1 is the more luminous one in the $G$ band (cf. Section 7.6.5).

Local optimisation

Once global optimisation has hopefully identified a set of parameters that place the cost ($\chi^{2}$) function inside the convex region of the global minimum, the final step is to perform a local optimisation in order to pinpoint the exact location of the global minimum. This is accomplished using the Levenberg-Marquardt algorithm which yields the parameter values that best match the observed LC, i.e., that minimise $F2$ (cf. Section 7.2.3). This is the most CPU-intensive part of the EB processing because Levenberg-Marquardt must call the EB simulator every time it needs to evaluate the cost function.

The following is a list of parameters that are fitted during local optimisation:

• Fill factors $s_{1}$ and $s_{2}$: The initial values are as provided by the global optimiser. During local optimisation, the fill factors are varied in a way that preserves the configuration determined by the global optimiser (detached, semi-detached or overcontact). In the case of an overcontact system, only one fill factor is fitted since the other one is determined by the common envelope constraint (cf. Section 7.6.3).

• Eccentricity $e$ and argument of periastron $\omega$: The initial values for both parameters are as provided by the global optimiser. If the global optimiser determined that the system is eccentric but the local optimiser in the end returns a value of $e<10^{-5}$ then $e$ and $\omega$ are both set to zero and the system is considered circular. However, if the global optimiser determined that the system is circular, i.e., that $e=\omega=0$, then the system is not allowed to become eccentric during local optimisation. See also Section 7.6.6.

• Inclination $i$: The initial value is as provided by the global optimiser. See also Section 7.6.6.

• Effective temperature $T_{\mathrm{eff},1}$: The initial value is as provided by the global optimiser. Since BP and RP data are not used in the fitting process (cf. Section 7.6.2), it is not possible to determine individual effective temperatures of the two components. However, a single-band LC (such as in $G$) can still constrain their luminosity ratio. Since the EB simulator cannot use luminosity ratios but requires temperatures as input, the strategy that was followed was to have the local optimiser fit $T_{\mathrm{eff},1}$ while keeping $T_{\mathrm{eff},2}$ fixed. Therefore, the individual temperatures at the end of the optimisation are not relevant or even correct, but the the luminosity ratio that they yield is.

Since in this case the temperatures are not physically relevant by themselves but only inasmuch as they affect the luminosity ratio, the EB simulator is invoked with a blackbody atmosphere. Indeed, there is nothing to be gained in this case from the added complexity of a more realistic model atmosphere.

It should be noted that the best-fit LC obtained from global optimisation is characterised by the “real” temperatures that were used in the simulated components. However, those temperatures can be arbitrarily rescaled in a way that still yields the same luminosity ratio (other physical considerations about the system notwithstanding), which is the only quantity that is actually constrained by the data.

• Time $t_{0}$ of eclipse of star 1 by star 2: The initial value is as provided by the global optimiser.

The following is a list of parameters that are required by the EB simulator during local optimisation but are neither fitted nor published as part of the solution:

• Orbital period $P$: The orbital/photometric period is provided by variability processing (cf. Section 7.6.2). See also Section 7.6.6.

• Mass ratio $q$ $(M_{2}/M_{1})$: The mass ratio is as provided by the global optimiser. It does not affect the LC when both fill factors are small, but it starts doing so the more the system approaches a semi-detached or contact configuration (i.e., as one or both fill factors become large enough). For Gaia DR3, it was assumed that the mass ratio determined in global optimisation is close enough to its true value that it can be taken as fixed during local optimisation. This is clearly an approximation, so the mass ratio is not among the published parameters. A better approach is being investigated for future data releases.

• Asynchronous rotation ratios $f_{1}$ and $f_{2}$: They are as provided by the global optimiser. Note, however, that when fitting an eccentric system, their value is recalculated via Equation 7.44 every time the simulator is called.

• Effective temperature $T_{\mathrm{eff},2}$: The initial value is as provided by the global optimiser. See description of $T_{\mathrm{eff},1}$ above for additional details.