© The Authors 2024

1. Introduction

The data of the space mission Gaia (Gaia Collaboration 2016) have stimulated investigations into open clusters (hereafter OCs) in the Galaxy. One of the areas of research that has significantly benefited is related to the binarity of OCs. While a lack of observational data of binary star clusters in our Galaxy in comparison with the Magellanic Clouds has been discussed previously (e.g. Vázquez et al. 2010; de La Fuente Marcos & de La Fuente Marcos 2009), recent investigations have provided a number of examples of pairs of OCs in our Galaxy that suggest a common origin (see, e.g. Qin et al. 2023; Angelo et al. 2022; Ye et al. 2022; Casado 2022a; Bisht et al. 2021; Zhong et al. 2019). Some clusters that were previously considered single have shown a binary signature upon detailed examination (Alejo et al. 2020; Dias et al. 2018), including clusters with signs of merging by tidal capture (Camargo 2021), according to Gaia data.

There is observational evidence that indicates binary star clusters are usually coeval and form as a result of hierarchical clustering in parental molecular clouds (Casado 2022b; Angelo et al. 2022). These pairs may or may not be bound and physically interacting. Notably, Angelo et al. (2022) found four bound and two interacting but unbound pairs of star clusters. Alternatively, coeval clusters may be located relatively close to each other due to being born in the same star forming complex but may never interact (Casado 2022a,b; Casado & Hendy 2023), and eventually they may become disaggregated while maintaining similar Galactic orbits (Coronado et al. 2022). On the other hand, pairs of closely located OCs with significantly different ages, which represent bona fide events of an occasional encounter, have been found (Camargo 2021; Piatti & Malhan 2022; Olivares et al. 2023). Dynamical simulation of such encounters has so far been limited to the publicly available orbit integration code galpy (Bovy 2015). The dynamical simulations of the evolution of initially bound primordial OCs in the Galactic tidal field by de la Fuente Marcos & de la Fuente Marcos (2010) have demonstrated that very close pairs, if formed, are merged in a short timescale, and in other pairs, either the smaller component is disrupted or a pair is decoupled in the background tidal field due to mass loss, making open star cluster binarity a transient phenomenon. Priyatikanto et al. (2016) used N-body simulations to explore the evolution of binary open star clusters in the Galactic tidal field with different initial orbital orientations with respect to the Galactic plane, which proved to be another important factor influencing the dynamical fate of such clusters. Darma et al. (2021) carried out N-body simulations to investigate the formation and evolution of binary star clusters in the Milky Way and in the Large Magellanic Cloud and found that binary open star clusters can form from stellar aggregates with a variety of initial conditions.

In Kovaleva et al. (2020), hereafter Paper I, we discussed possible scenarios of origin and dynamical evolution in a simplified approach and constrained parameters, within the limits of observational values, that suggest a common origin of the open star clusters Collinder 135 (hereafter Cr 135) and UBC 7 (see also discussion in Pang et al. 2021). As one of the results of that work, we identified, using backward orbital integration, a presumable grouping of initial states in 6D space that corresponds to a common origin of these clusters. Our previous work was limited only by the investigation of motion of OC dynamical centres (centre of masses). So, we completely neglected the OC internal dynamical and stellar evolution, including the mass loss of the clusters. By overcoming such restrictions in our current investigation, we make our previous results more robust and reliable.

In this article, our purpose is to reproduce the dynamical evolution of the system of Cr 135 and UBC 7 using star-by-star simulations. During our investigation, we especially paid attention to the complex dynamical co-evolution and cross-penetration of the clusters’ stellar populations. This behaviour of the modelled stellar systems made our attempts of system modelling quite challenging.

In comparison with Paper I, we have updated the numerical implementation of the φ−GPU with a new stellar evolution block (Banerjee et al. 2020; Kamlah et al. 2022) combined with our high-order dynamical N-body code. This new implementation allowed us to predict with a good degree of confidence the evolved clusters’ internal stellar structure and cumulative mass distribution.

After the publication of Paper I, the Gaia DR3 catalogue (Gaia Collaboration 2021, 2023) was published. With respect to the previously used Gaia DR2 catalogue (Gaia Collaboration 2018b), Gaia DR3 provides enhanced astrometric and photometric data for an enlarged sample of the stars as well as a significantly larger number of sources with mean radial velocities. Some systematic effects have been resolved or mitigated regarding Gaia DR2, which is also beneficial for the study of the population of OCs. These improvements allowed us to obtain refined characteristics of the system composed of the two clusters based on the deep census of the probable members of this system.

Our goal is to reproduce the present-day state of OCs and compare it with the most up-to-date information about their stellar population from Gaia DR3 data. Such an ambitious investigation has only been possible since last few years thanks to the Gaia’s unprecedented position and proper motion measurements. As a possible physical pair of OCs, Cr 135 and UBC 7 give us a unique laboratory to also study the common processes of the ongoing star formation on larger scales in our Milky Way disc.

The paper is organised as follows. In Sect. 2, we describe dynamical N-body modelling of evolution of the pair of clusters. In Sect. 3, we discuss the data filtering process that resulted in the following view onto the system of the two OCs, Cr 135 and UBC 7, with Gaia DR3 data. Section 4 is focused on characterisation of the system based on the obtained characteristics of its members. In Sect. 5, we describe results obtained in our N-body simulations in comparison with the observational data presented in Sects. 3 and 4. In Sect. 5, we also discuss kinematics of the surrounding extended halo. In Sect. 6, we sum up the main results of the work.

2. Cr 135 and UBC 7: Star-by-star N-body modelling

2.1. Initial conditions and integration procedure

In Paper I we found that it is highly probable that Cr 135 and UBC 7 form a physical pair. This conclusion is based on backward-in-time numerical simulations and 6D observations from the Gaia DR2 catalogue. The resulting initial positions and velocities for cluster centres of mass (initial separation and velocity: 10 pc and 0.95 km s⁻¹) showed that the probability of a mutual formation scenario is about 98% (see the numerical model #[53,61], and Fig. 2 of Paper I). In order to estimate the impact of Gaia DR3 data on new simulations, we used newly determined present-day cluster observed parameters determined as described in detail in Sect. 3. We found that there is good agreement in the newly determined Galactic coordinates, proper motions, and radial velocities (see Table 1) with the Gaia DR2 parameters derived in Paper I.

To evaluate the possible distinction (DR2 vs. DR3) in the distances between the two OC centres, we carried out a new lookback integration using the new DR3 initial data from Table 1 and compared it to our earlier selected best numerical model #[53,61] (see Fig. 1). As can be seen in the figure, the agreement is good, and the effect of the different data releases is small.

Fig. 1. Evolution of the distance between the centres of Cr 135 and UBC 7 for the DR2 (red) and DR3 (blue) initial conditions.

For the OCs’ dynamical orbital integration including the stellar evolution, we used the high-order parallel N-body code φ−GPU¹, which is based on the fourth-order Hermite integration scheme with individual hierarchical block time steps (Berczik et al. 2011, 2013). The current version of the φ−GPU code uses a native GPU support and a direct code access to the GPU using the NVIDIA native CUDA library. The present code is well tested and was used to obtain important results in our earlier OC simulations Just et al. (2009), Shukirgaliyev et al. (2017, 2018, 2021).

Because of the relatively short evolutionary timescale of our system, we used the so-called fixed external potential with parameters independent of time. Similar to Kharchenko et al. (2009), Ishchenko et al. (2023), we used the three-component (bulge-disc-halo) axisymmetric Plummer-Kuzmin model (Miyamoto & Nagai 1975).

The current code also takes into account up-to-date stellar evolution models. The most important updates were made to components of stellar evolution and include updated metallicity dependent stellar winds; updated metallicity dependent core-collapse supernovae, along with a new fallback prescription and their remnant masses; updated electron-capture supernovae accretion induced collapse and merger-induced collapse; remnant masses and natal kicks; and black hole natal spins, among other updates. (For more details about the new stellar evolution library, see Kamlah et al. 2022).

To carry out the star-by-star N-body modelling together with the stellar evolution for both OCs, we assumed the clusters to be in initial dynamical equilibrium. For the star particles, we used the King model distribution function King (1966) with three free parameters: the initial total mass, the initial half-mass radius r_hm, and the dimensionless central potential W₀. Using the Kroupa (Kroupa 2001) initial mass function with the lower and upper mass limits equal to 0.1–10 M_⊙, and the classical solar value for the initial stellar metallicity, Z = Z_⊙ = 0.02 (Grevesse & Sauval 1998), for both clusters. For the centres of the clusters, we used the initial coordinates and velocities at the 50 Myr lookback time integration (see Table 2).

We tested more than 50 runs with various combinations of M, r_hm, and W₀ separately for each cluster using the full N-body model calculations up to the present time. For the initial positions and orbital velocities in Cartesian Galactic coordinates of the star clusters (SC) centre at the −50 Myr lookback in time, we used the values obtained during our own set of lookback integration of the SC centre of mass. The initial masses, radii, and King profiles were found during the fitting test runs, with the aim of matching the theoretical and observed parameters. The basic initial parameters of the model are given in Table 2 separately for each SC.

Figure 2 shows the observed stars of the OCs (light green for Cr 135 and light blue for UBC 7) overplotted with the present-day positions of all particles from our basic model according to the initial data in Table 2 (dark green for Cr 135 and dark blue for UBC 7). The mutual orientation of the OCs is correct, although both sets of particles (representing Cr 135 and UBC 7) are slightly misaligned with their prototype clusters (see Table 3). In Table 3, we present the positions and velocities for the current centres of mass for both OCs. The identifier Sim indicates the cluster’s centre of mass as obtained from our basic model, and Obs indicates the centre of mass as obtained from observations (see Sect. 3). The small disagreement between the 3D coordinates and velocities inside the table most probably arises from the different mutual interactions of the clusters in the backward integration of the orbit compared to the forward integration of the full OCs. The backward integration was performed with the present-day masses of the OCs. The forward integration started with the larger initial cluster masses and includes mass loss by stellar evolution and the tidal field. The present-day shift of the cluster centres demonstrates that Cr 135 and UBC 7 interact gravitationally with each other.

Fig. 2. Comparison of observed positions of probable cluster members and of model stars of basic N-body simulations. The light green and light blue dots are the members of Cr 135 and UBC 7, respectively. The coordinate frame was obtained by parallel translation of the standard Galactic coordinate system.

2.2. Stellar evolution and lifetime

In Fig. 3, we present the evolution of the tidal radius and bound particle numbers according to our basic N-body model. The radii of both clusters show a slight increase during the first several tens of million years, after which the radius of Cr 135 starts to sharply decrease, and after approximately 110 million years of evolution, the core falls apart (technically, the number of bound particles decreases to zero). At this time, the separation between the two OCs is about 50 pc. UBC 7 has a more gradual decrease in tidal radius in contrast to Cr 135, and it falls apart approximately 50 million years after the dissolution of Cr 135. The differences between the dynamical behaviour of our clusters can be explained as being a result of the concentration parameters W₀ = 3 and W₀ = 11 of the different initial King profiles. Thus, for Cr 135 the profile indicates a looser and more unstable core. In UBC 7, the core is more concentrated, and therefore, despite the smaller number of stars, it demonstrates a greater ability to survive.

Fig. 3. Evolution of the tidal radius (upper panel) and bound particle number (bottom panel). The green line represents Cr 135, and the blue line indicates UBC 7.

In Fig. 4, we show the orbital evolution up to ∼150 Myr from the formation of the OCs. As can be seen from Fig. 4, during the time of evolution, the clusters follow mildly eccentric clockwise orbits around the Galactic centre, as seen from the Galactic north pole, and simultaneously rotate around each other. The orbits have considerable vertical oscillations in the z-direction between ±80 pc.

Fig. 4. Orbital evolution of Cr 135 and UBC 7. The colour-coding represents the time of integration. The filled and open black circles indicate the present positions of the clusters. Th black crosses indicate the positions where the clusters dissolve.

In Fig. 5, we present the Hertzsprung-Russell diagram (HRD) for three moments of evolution (50, 100, and 150 Myr in forward integration) for each OC starting from 0 Myr. Due to the low number of available stars and the short period of the OCs’ lifetime, we mainly have main-sequence stars. When analysing the type of stars for the present day (which have an age of 50 Myr according to the N-body simulation), ∼84% of them are low-mass stars (m < 0.8 M_⊙), and only ∼16% of them are higher-mass main-sequence stars for both OCs. At the end of the forward integration, low-mass stars were still dominant, and there were only a few carbon/oxygen and oxygen/neon white dwarf stars.

Fig. 5. Hertzsprung-Russell diagrams of Cr 135 and UBC 7 (including all unbound stars). Red dots correspond to an age of 0 Myr, light blue to 50 Myr, green to 100 Myr, and black to 150 Myr. The crosses and the stars show two types of white dwarfs from the N-body simulation.

3. Cluster membership with Gaia DR3 data

Our work in Paper I is based on Gaia DR2 data and aimed at selecting the most probable members of Cr 135 and UBC 7 in order to obtain their global space and kinematic parameters. We incidentally discovered an extended corona of stars with an almost constant projected density inside a 20-pc radius around the centre of mass of both clusters. In the present investigation, our aim is to use the wealth of the Gaia DR3 data to explore stellar membership and the structure of the whole system in detail.

Recent investigations have proved the existence of extended features connected with groups of OC (see, for example, Beccari et al. 2020; Zhong et al. 2019; Pang et al. 2020). On the other hand, the discussed Cr 135 and UBC 7 are situated in the Vela-Puppis star-formation region (Galactic longitude 245° ≲l ≲ 265° and latitude −15° ≲b ≲ −5°), which hosts several young OCs with ages between approximately 15 and 50 Myr (Cantat-Gaudin et al. 2019b,a; Beccari et al. 2020). Of these clusters, Cr 140, NGC 2451B, and NGC 2547 are situated at about 100 pc from the investigated pair. The nearest neighbour to Cr 135 and UBC 7 is the evolved open cluster Alessi 3 (625 Myr, Hunt & Reffert 2023). Thus, in order to reveal the size and characteristics of the coronae and check whether it is due to contamination from neighbouring stellar groups, we considered a very wide region surrounding Cr 135 and UBC 7. This required refining of the data selection principles used in Paper I.

We performed an intensive search for members related to OCs Cr 135 and UBC 7 among all Gaia DR3 sources satisfying the following conditions:

$$\begin{array}{c}\hfill \begin{array}{c}2.5<\varpi [\mathrm{mas}]<5,\hfill \\ 220°<l<280°,\hfill \\ -25°<b<1°,\hfill \\ 1.2<{\mu }_{\alpha }^{\ast }[\mathrm{mas}{\mathrm{yr}}^{-1}]<11.8,\mathrm{and}\hfill \\ -15<{\mu }_{\delta }[\mathrm{mas}{\mathrm{yr}}^{-1}]<-5.\hfill \end{array}\end{array}$$

The linear size of a region was determined based on previous research in order to include regions that are most probably free of members of central clusters. These external regions are useful for the estimation of random background contamination from the procedure. We applied parallax zero-point correction (Lindegren et al. 2021) with the algorithm provided by Gaia².

Figure 6 presents the primary selection containing 116 975 sources. Their colour-coded parallax ranges from yellow, for the most distant stars of the sample at about 400 pc from the observer, to dark violet, for foreground stars at 250 pc from the Sun. Even with this rough selection, some of the OCs of a region became noticeable, Cr 135 and UBC 7 being the nearest two. At a distance of about 280–300 pc from the Sun where they are situated, 1 deg ≈ 5 pc, so the linear size of the rectangle in the frontal plane is about 300 pc in l and about 75 pc in b (Fig. 6, left panel). One can see an enhanced density of sources towards the Galactic plane in b and towards the Vela Molecular Ridge and the Galactic centre in l. As we applied no filter to the quality of astrometric and photometric solutions, part of the faint stars may belong to the background and may reflect the stellar density of distant regions (Lindegren et al. 2018).

Fig. 6. A selection of Gaia DR3 data in the broad vicinity of the Cr 135 and UBC 7. Cr 135, UBC 7, and their closest neighbouring clusters according to Cantat-Gaudin et al. (2020), left panel. Refined sample of 60% of the members of Cr 135 and UBC 7 over the same sky area, right panel. The colour-coding indicates the MP according to Eq. (1). The circles in both panels represent the half-number radii of the clusters according to Cantat-Gaudin et al. (2020). See discussion in the text.

We kept as a basis the procedure of probabilistic selection of cluster members described by Kharchenko et al. (2012). Each star of the considered dataset was assigned a probability P of sharing age, space, and kinematics of a pair of clusters, namely Cr 135 and UBC 7:

$$\begin{array}{c}\hfill P=min(P_{\mathrm{kin}},P_{\varpi },P_{\mathrm{ph}}).\end{array}$$(1)

where P_kin, P_ϖ, P_ph are the probabilities based on kinematics, parallax, and photometric data. The procedure of calculating the probabilities is described in detail in Appendix A.

After each source in our dataset was assigned a membership probability (MP) according to Eq. (1), we selected stars with P > 0.6 (so-called 60%-members) for further consideration. In Fig. 6 (right panel), one can observe two dense central groups (i.e. the cores of Cr 135 and UBC 7) and an extended group of stars with a high MP resembling a halo surrounding the cores. Also present are sparse, chaotically distributed sources of lower probability that dominate the periphery, which presumably are contaminating foreground and background stars, and an elongated dense structure between 250° < l < 263° and −11° < b < −7° formed by members with an intermediate MP. After investigation, we found that this elongated stellar group has a systematically different space motion and does not belong to the system we are focused on. This group may in part be associated with a group of stars identified as a loosely spread cluster [BBJ2018] 6 in Beccari et al. (2018), but due to its spatial scale, it might be referred to as a stream instead. For the purposes of the present work, we identified and removed from consideration sources with MPs to the [BBJ2018] 6 group larger than MPs to the group of Cr 135 and UBC 7. The remaining dataset represented in Fig. 7 consists of 1714 stars. One may observe that the probable members of the group composed of the two clusters still form a distinct asymmetrical halo around the location of the dense clusters.

Fig. 7. Components of the OC system. The probable members of Cr 135 are shown with red dots, probable members of UBC 7 are indicated with blue dots, and extended halo members are shown with yellow dots. The upper panel shows the lb Galactic coordinates sky plane; the bottom panel shows the ϖ,  b plane.

Random contamination can be estimated based on the number of probable members in the outer parts of the field, which are approximately at l > 270°, l < 235°, b > −4°, b < −20°. For this purpose, as well as to set a distinction between the halo and the background, we relied on the density distribution of the stars in the tangential plane. For each source, we estimated its weight over the background W_i as the residual:

$$\begin{array}{c}\hfill W_i=\rho (l_i,b_i)-{\rho }_{\mathrm{r}}(l_i,b_i),\end{array}$$(2)

where ρ(l, b) is the stellar density distribution of the refined dataset convolved by the 2D Gaussian kernel of the width of 2° and ρ_r(l, b) is the similarly convolved random flat distribution of the same number of stars over the same field. We considered the 293 sources with W_i ≤ 0 to be background contamination stars, and we treated the 1421 stars with a positive W_i as those composing the clusters and the halo. From hereon, we refer to this dataset as Sample I (available via CDS).

Our purpose was to recover the most probable distribution of sources in Sample I between three substructures: Cr 135, UBC 7, and their extended halo. In space, these groups appear to be separated, but there is no distinct border between them. Given that the separation between the centres of the clusters is just 22 pc, peripheral stars of the clusters are mixed and may only be conditionally attributed to one of them. Due to the similarity of their mean kinematic characteristics, we attributed each source to one of these components based solely on its maximum likelihood location in space with respect to the clusters as follows.

Using the mean position of the reliable members from Table 1, we suggested for the clusters a Gaussian 3D PDF $P_k^{\text{sc}}$ over a Cartesian distance to the centre with coordinates based on Gaia DR3 ($l_k^0,b_k^0,{\varpi }_k^0$) from Table 1, where k = 1 refers to Cr 135 and k = 2 is for UBC 7:

$$\begin{array}{c}\hfill \mathrm{\Pi }(d_k)=exp\left[\frac{{\left(d_k^i\right)}^2}{2{\sigma }_{d^k}^2}\right],\end{array}$$(3)

where σ_d^k is based on the distribution of the reliable cluster members by Gaia DR3 (Table 1). We considered the distance d between the centre of an OC ($l_k^0,b_k^0,{\varpi }_k^0$) and the i^th source of Sample I (lⁱ, bⁱ, ϖⁱ) to be a function of ${\varpi }_m^i$ with a Gaussian PDF:

$$\begin{array}{c}\hfill \mathrm{\Pi }({\varpi }^i)=exp[-\frac{{({\varpi }^i-{\varpi }_m^i)}^2}{2{\sigma }_{{\varpi }^i}^2}],\end{array}$$(4)

where the standard deviation ${\sigma }_{\varpi }^i\equiv {ϵ}_{\varpi }^i$. We simulated numerically the resulting PDF for the distance of each source of Sample I to the cluster centre $\mathrm{\Pi }(d_k,{\varpi }^i)=\mathrm{\Pi }(d_k)\cdot \mathrm{\Pi }({\varpi }_k^i)$, and we compared the maximum values of PDFs for k = 1, 2:

$\max (\mathrm{\Pi }(d_1,{\varpi }^i))\ge \max (\mathrm{\Pi }(d_2,{\varpi }^i))\wedge d_1^i[\max (\mathrm{\Pi }(d_1,{\varpi }^i))]\le 15$ pc, then

   the source i∈ Cr 135,

else

   if $\max (\mathrm{\Pi }(d_2,{\varpi }^i))\ge \max (\mathrm{\Pi }(d_1,{\varpi }^i))\wedge d_2^i(\max (\mathrm{\Pi }(d_2,{\varpi }^i)))\le 15$ pc, then

     the source i∈ UBC 7,

   else

     the source i∈ halo.

   end if

end if

The border between the cluster and halo members was selected to be equal to 15 pc, tentatively, based on the results of the analysis of the spatial distribution of Sample I. As a result, we related 294 stars to Cr 135, 243 stars to UBC 7, and 884 stars to the extended halo (see Fig. 7).

4. Estimation of the parameters of the clusters

The system includes two quasi-ellipsoidal dense stellar groups that we refer to as clusters, which have similar but distinct mean space motions. They are surrounded by an extended stellar halo of irregular form and have an approximately constant stellar density in projection onto the lb-plane (Fig. 7, upper panel). This structure extends for more than 80 pc towards the Galactic anti-centre and half as much in the opposite direction. The dimensions towards and from the Galactic plane are comparable with the size of the clusters. In space, the structure is aligned along the direction joining the centres of Cr 135 and UBC 7, becoming closer to the Galactic plane the farther away it is from the Sun (see Fig. 7, bottom panel).

The similarity of the populations of two clusters and of the population of the extended halo is illustrated by their distribution over the colour-magnitude diagram in Fig. 8. In the figure, apparent magnitudes and reddened colours have been reduced to the absolute scale using individual parallaxes of the stars and mean reddening of the central probable members as in Paper I. The Green line represents the Padova isochrone for 50 Myr for single stars, and the magenta line is for the same isochrone for the unresolved binary stars with identical components. The red line is for the zero-age main sequence built in Paper I as a hot envelope of the related set of Padova isochrones of different ages. The probable members of the components of the system are represented by the same colours as in Fig. 7. For reference, grey dots represent each tenth background source of the same spatial region (as in Fig. 6). One may notice not only the similar age of the OCs and the halo but also the similar density distribution of probable members along the isochrone, which suggests a common origin.

Fig. 8. Components of the OC system at the colour-magnitude diagram. Probable members of Cr 135 are shown with red dots, probable members of UBC 7 are indicated with blue dots (left panel), and probable members of the extended halo are shown with yellow dots (right panel). Grey dots are for background stars. The green line represents the isochrone for 50 Myr for single stars, and the magenta line is for the same isochrone for unresolved binary stars with identical components. The red line is for the zero-age main sequence.

The characteristics of the components of a system composed of the two OCs and their halo based on Sample I data selection are listed in Table 4. This includes the number of probable members as defined in Sect. 3, N, and the half-number radii of the clusters and halo, R_HN.

One may expect that the halo contains significant contamination and a random background. Its numerical estimate based on the surface density of a random background in the neighbouring regions led to 0.25 ± 0.05 contaminating stars per square degree. This decreased the expected actual population of the halo by 60–120 stars, taking into account the vagueness of its borders.

4.1. Counted masses

The estimate of the present-day mass of the system may be assessed based on masses of visible probable members as ≈900 M_⊙ (see Table 4, M^Σ). More than half of the mass lies in the halo. There are a few probable causes of it being underestimated, such as hidden mass in low-mass stars beyond the photometric limit, members with an MP less than 60%, and unresolved binarity. As well, there is likely random contamination that causes present-day mass overestimation. The effect of random contamination is mostly important for the halo due to its large extent. We provide a numerical estimate for the expected contamination based on the surface density of the random background in the neighbouring regions as obtained in Sect. 4. Filtering out the expected background reduced the mass of the halo by 40–70 M_⊙. For the clusters, the halo itself should be considered a contaminating background; the estimate of the mean density for it is 3.6 ± 0.2 stars per square degree. This provides an upper estimate for a contaminating mass of ≈27 M_⊙ for Cr 135 and ≈17 M_⊙ for UBC 7.

A possible source of mass underestimation lies within the unresolved stellar multiplicity and, primarily, binarity (see discussion in Borodina et al. 2019, 2021). The effect depends on the fraction of binary to multiple stars and their distribution over the mass fraction q = m_B/m_A, where m_A and m_B are respectively masses of the most massive and least massive components of a binary star. We estimated the fraction of binary stars α = N_bin/N_tot with mass ratio of components q > 0.4 for the system directly from the colour-magnitude diagram (Borodina et al. 2019; Borodina & Kovaleva 2020) as α = 23%±12%. The result is consistent with the binary fraction estimates carried out by Pang et al. (2023). The estimated fraction of binary stars leads to the increase of the expected mass of the system by a factor of 1.1.

The complete mass of all stellar members of the system includes low-mass stars beyond the photometric completeness limit. We estimated the total hidden mass, including the unresolved components of binaries, based on the suggestion of a Kroupa-like IMF (Kroupa 2001) of the system, with normalisation to the number of stars with a visible G magnitude phot_g_mean_mag_corr < 17.5. This provided an estimate of the total mass, including unseen stars with masses down to 0.09 M_⊙, that is, 1.3–1.4 times larger than the straightforward sum of masses of observed sources. These suggestions provide a scale of reliability of counted mass and resulted in the asymmetrical uncertainties indicated in Table 4.

In Table 4, M^Σ is the apparent counted masses of the 60%-probability members. The indicated confidence limits in the table represent an extreme possible increase and decrease of this value due to taking into account hidden members or filtering out a maximum possible random contamination. The mean mass of a star is shown by $\overline{M_i}$. The mean mass of the probable members of Cr 135 $\overline{M_i}$ exceeds the mean mass of the probable members of UBC 7 and that of the halo (which are similar; see Table 4). This may be a result of mass segregation.

While the halo as a whole is more massive than both clusters together, its density is much lower, below 0.001 M_⊙ pc⁻³. The mean stellar density within 15 pc from the centre of Cr 135 and UBC 7 is 20–30 times larger, and in the centre of the cores of the clusters, it is up to around 3.5 M_⊙ pc⁻³.

4.2. Tidal parameters of substructures of the Cr 135 and UBC 7

The spatial structure of an OC can be standardised with the help of a well-known empirical model proposed by King (1962). The model parameterises the radial surface density profile of a cluster in terms of its core and tidal (total) radii r_c and r_t. In the integrated form, the profile can be represented with the following equation:

$$\begin{array}{cc}\hfill n(r)=\pi r_{\mathrm{c}}^{2}k& \{ln[1+{(r/r_{\mathrm{c}})}^2]-4\frac{{[1+{(r/r_{\mathrm{c}})}^2]}^{1/2}-1}{{[1+{(r_{\mathrm{t}}/r_{\mathrm{c}})}^2]}^{1/2}}\hfill \\ \hfill & +\frac{{(r/r_{\mathrm{c}})}^2}{1+{(r_{\mathrm{t}}/r_{\mathrm{c}})}^2}\},\hfill \end{array}$$(5)

where n(r) is the number of stars projected onto the ideal (tangential) plane within a circle with radius r and k is a normalising coefficient.

The observed profiles were counted based on the sources of Sample I with brightness limitation phot_g_mean_mag_corr < 18.0. In order to avoid mutual contamination, we evaluated the residual background level produced by the neighbouring population via estimation of the average surface density produced by the sample stars in the neighbouring circular area outside the visual limit of the respective cluster. The fit and determination of the parameters was carried out using the MPFIT procedure from the Markwardt (2009) IDL library.

Figure 9 shows the empirical profile fits for our working samples. Figure 10 illustrates how reasonably computed tidal parameters outline the general structure of our pair of star clusters. As one can see, the core radii clearly mark off the densest parts of the star clusters. In 3D, both clusters are touching at their tidal limits. We regard this as evidence of the physical relation of both clusters residing in the common potential well. In fact it is clear that the bodies that are currently recognised as separate OCs themselves are rather double cluster cores immersed in the common halo.

Fig. 9. Fits of the King profile according to Eq. (5) for the observed radial density profiles, based on the confident members of the clusters. The curves present distributions including the residual background stars (black dots on the yellow lines), the distributions corrected for the background (circles with statistical error bars), and the fitted King profiles (solid red curves). The vertical dashed lines show the core and tidal radii (rc and rt). The colour bands show errors in the radii.

Fig. 10. Open clusters Cr 135 and UBC 7 in the tangential plane contacting the celestial sphere at the common centre of the mass of the clusters. The red dots show probable members of Cr 135, and blue dots show probable members of UBC 7. The circles show the core and tidal radii. Black dots show the centres of the clusters. The “+” sign shows the centre of mass, which is assumed to be at 292 pc from the Sun. Linear coordinates ξ and η are parallel to equatorial ones.

In Table 4 we show the derived profile parameters. The fit procedure works well in our case, and the profile correction for residual outer background plays a critical role for the correct determination of the tidal radii.

In order to understand whether one can observe a footprint of cluster gravitational interaction in the observed structure of the pair of star clusters, we tried to compare tidal parameters for the component profiles counted in sectors opened at different directions. As we found no systematic difference is observed, the derived tidal parameters vary within the fit errors only.

Knowledge regarding cluster structure tidal parameters can be used for independent mass estimates of an OC. To determine the tidal mass m_t, we followed King (1962), Ernst et al. (2010) and Just et al. (2023) and used a condition for the balance of gravitational forces between the Galaxy and the OC:

$$\begin{array}{c}\hfill m_{\mathrm{t}}=\frac{r_{\mathrm{J}}^3}{G}(4-{\beta }^2){\mathrm{\Omega }}^2,\end{array}$$(6)

where r_J is the Jacobi radius (distance from the cluster centre to Lagrange points L₁, or L₂, where its own gravity is equal to the Galaxy field), Ω is the angular velocity of rotation of the Galaxy at the cluster’s galactocentric distance R, β is the ratio of epicyclic κ and rotational Ω frequencies, and G is the gravitational constant.

The relation between Jacoby and tidal radii for realistic clusters in the case of the three-component Plummer-Kuzmin model of Galactic potential (Miyamoto & Nagai 1975) was studied using N-body calculations by Ernst et al. (2010). They showed that the ratio r_t/r_J depends on the cluster position on the sky (mainly on Galactic latitude) and on their age. Its typical value at low latitudes of |b|< 20° computed by Ernst et al. (2010) for age of cluster model of 0.62 Gyr is lower than r_t/r_J ≈ 1.19. For the location of our clusters at R ≈ 8.29 kpc, the rotation parameters computed from the Plummer-Kuzmin model are equal to β ≈ 1.45, Ω ≈ 28.20 km s⁻¹ kpc⁻¹. The respective values of the tidal masses are shown in Table 4, and the formal uncertainty was estimated via propagation of the error in r_t. Since our clusters are much younger than the lower age limit of the Ernst et al. (2010) model and since we observed that the ratio tends to decrease at younger ages, we believe that the found values manifest the upper limits of the respective tidal masses.

The mass estimate obtained for UBC 7 via correction of counted mass for hidden mass (≈204 M_⊙; see Table 4) agrees rather well with the tidal mass (270 ± 140 M_⊙). The tidal mass for Cr 135 (770 ± 340 M_⊙) exceeds the estimate obtained from the counted mass even after correcting for the expected hidden mass (≈305 M_⊙) and taking into account the large uncertainty. This may be related to difficulties with the treatment of a populated halo. In addition, as noted previously, the found values probably manifest the upper limits of respective tidal masses.

One may also compare present-day tidal radii with the results of the dynamical evolution modelling (Fig. 3 for the present time T = 50 Myr). The model predicts very similar tidal radii for the two OCs, about 12 pc for Cr 135, which is slightly less than for UBC 7. Taking into account the errors, the results match the tidal radii obtained from observations quite satisfactorily, especially for UBC 7. For Cr 135, the estimated tidal radius (r_t = 15.6 ± 2.6 pc; Table 4) is larger than what the dynamical evolution model predicts.

5. Discussion

5.1. Cumulative star number and proper motion distributions inside the modelled clusters

In this section, we present a comparison of the present-day radial cumulative star count obtained from our basic model with initial conditions presented in Table 2 with the current observational data (see Sect. 4.2, Fig. 9 with Table 4). In Fig. 11, the total cumulative number distribution (CND) of cluster stars excluding the stellar background but including the observational errors in number count are presented as a grey zone. For Cr 135 and UBC 7, we used the fit range within 0 < r < 20 pc. We set such a large distance limit in order to exceed the OCs’ current Jacobi radii (see Table 4). We present our basic model as a thick black line in these plots, which was obtained directly from N-body simulations. Because of the magnitude limit and incompleteness of the faint sources in the observational data due to the Gaia satellite specifications, we selected from the numerical models only the stars that are within the specific observed stellar mass range: from 0.28 to 4.0 M_⊙. In this case, our model (black line) has been adopted as an observed CND of stars for both clusters.

Fig. 11. Cumulative number distribution of stars compared against the projected radius for Cr 135 (left panels) and for UBC 7 (right panels). The coloured lines represent different N-body randomisations. The upper panel presents the IMF sampling, and the bottom panels show the PV sampling. The shaded grey regions show the observational data with uncertainties.

After the physical parameters of the model clusters were fixed, we carried out an additional set of 2 × 9 = 18 random numerical realisations of the same physical model. In the first nine realisations, we used different samples of the initial mass function (IMF) of the clusters, while the initial positions and velocities of the stars were fixed. In the second set of nine realisations, we randomised the initial positions and velocities (PVs) of the stars, while the IMF distribution was fixed. This enabled us to estimate the possible uncertainties directly connected with the randomisation of the initial individual stellar masses and PVs phase-space distributions.

In Fig. 11, all the coloured lines are inside the grey area, especially in the central 8 pc, which is where 80% of the stars reside. However, our plots show that half of the coloured lines are on the edge of the grey area, around the 10–11 pc. We interpret this as an indication of a somewhat smaller concentration of the real initial OCs’ PV model distribution compared to our initially selected numerical values. Taking into account the OCs’ centres of mass uncertainties obtained from the 18 samples, we got the next standard deviation in the three Cartesian Galactic coordinates for Cr 135 (±2.70, ±2.50, and ±1.86 pc) and for UBC 7 (±1.85, ±1.41, and ±0.96 pc). As can be seen, we got slightly larger deviations for Cr 135 as compared to UBC 7.

In Fig. 12, we present the same CNDs but compare them against the residual proper motion magnitude inside the clusters. As can be seen from the plots, the proper motion velocity components of the stars are also quite consistent with the observed distributions and error bars. Together with the distance CNDs, these results show that our numerical modelling is able to reproduce the currently observed 6D phase-space distributions with a high level of confidence.

Fig. 12. Same as in Fig. 11 but compared against the residual magnitude of the proper motion.

5.2. Origin of the extended halo

In Sect. 3, we pointed out the presence of a halo surrounding both clusters with a mean density of 3.6 ± 0.2 stars per square degree. The halo stars are not in our current N-body simulations because we started with two King spheres in virial equilibrium that are accurate according to mutual attraction and the Galactic tidal field.

An explanation for the halo can be found in the cluster formation scenario with low star-formation efficiency (SFE; see Shukirgaliyev et al. 2021, and references therein). It suggests that most of the gas is expelled abruptly from the star-forming region due to stellar feedback from massive OB stars. A sudden loss of approximate equilibrium in the gas-stars system leads to so-called violent relaxation followed by star ejection.

To find evidence in favour of this scenario, we evaluated the net residual velocities of the stars with respect to the centre of the pair of star clusters in the line-of-sight and tangential directions. The obtained net line-of-sight velocities of the receding and approaching stars are shown in Fig. 13. In the figure, the mean residual radial and tangential velocities of the cluster members (red, blue, and green lines) are compared with low-SFE cluster values (grey line), which were obtained using data by Shukirgaliyev et al. (2021). In this paper, we used the Dehnen γ = 0 initial stellar profile with SFE = 0.05 after 50 Myr of evolution, which is presented for the first time. This choice of the model was determined by the bound mass of the numerical cluster, which should be comparable to the mass of our pair of star clusters. We calculated the runaway velocities of stars from the centre of mass of the pair of clusters with the projection on the sky plane. The mean value of the velocities is shown with the grey line in Fig. 13, and the shadow represents its standard deviation. The agreement between the net tangential observed velocity and the mean radial velocity in the simulations is excellent, supporting the scenario of low-SFE formation of the pair of star clusters.

Fig. 13. Outwards motion of the corona stars with respect to the pair of star clusters’ centre of mass. Smoothed magnitudes Vr of residual line-of-sight velocities of receding (approaching) stars from Gaia DR3 are shown by red (blue) lines. The mean residual tangential velocities in the picture plane Vt are shown by the green line. The grey line shows the mean tangential components of the runaway velocities obtained in N-body simulations for the Dehnen γ = 0 initial stellar profile with SFE = 0.05 after 50 Myr of evolution (Shukirgaliyev et al. 2021), where the shadow represents the uncertainty.

6. Conclusions

We performed star-by-star dynamical modelling of the evolution of the two OCs Collinder 135 and UBC 7 from their supposed initial state (mean coordinates and velocities) to their present-day state. The models were compared with the data obtained based on the performed census of membership of the two clusters with Gaia DR3 data. In particular, we used CNDs of the cluster members over the distance from their centres from observations to adjust the initial parameters of the clusters, such as the number of stars and their concentration parameter in the model. We selected these parameters in order to reproduce well the present-day state and distribution of the members of Cr 135 and UBC 7 in 6D space (coordinates, velocities).

Our research supports the potentially shared origin of the Cr 135 and UBC 7 OCs. By integrating the dynamical co-evolution of these clusters within an external Galactic potential, we obtained results that align closely with the internal stellar density and velocity distributions of the clusters. In our parameter-fitting procedure for the clusters, we focused solely on the resulting stellar CNDs (Fig. 11) at 50 Myr post cluster formation. The model’s resulting profiles for the residual proper motions and the Gaia DR3 population of the clusters show a strong agreement without the need for any additional fitting, Fig. 12.

The star-by-star dynamical evolution model shows the future of the stellar populations of Cr 135 and UBC 7 (including the existence of the stellar remnants, such as C/O and O/N white dwarfs). Due to the limited number of stars (only a few hundred) in both clusters, their 50 Myr-evolved stellar populations have a lack of neutron stars and black holes.

Based on our comprehensive dynamical modelling with the up-to-date stellar evolution prescription, we are able to predict the future ∼100 Myr evolution of both clusters with a good degree confidence. In our models, the UBC 7 cluster lives significantly longer compared to the Cr 135 cluster, Fig. 3. This is not surprising given the fact that, according to our numerical modelling, UBC 7 is denser (see Table 2) than Cr 135.

During our investigation, we found strong evidence of the presence of the common stellar halo in the Cr 135 and UBC 7 OC system, and we estimated its characteristics. The members of the halo demonstrate net residual velocities of stars with respect to the centre of the pair of star clusters that are compliant with the mean radial velocity of particles of N-body simulations of a low-SFE cluster stimulated by violent relaxation after gas loss.

Detailed common dynamical modelling of the OCs and the extended halo system was beyond of the scope of our current paper. Such modelling, however, will definitely be of interest and valuable for future investigations. The presence of such an extended halo can also have an impact on the clusters’ cumulative mass distributions and internal proper motions.

Acknowledgments

The authors thank the anonymous referee for a very constructive report and suggestions that helped significantly improve the quality of the manuscript. The work of MI, PB and OS was partially supported under the special program “Long-term program of support of the Ukrainian research teams at the PAS Polish Academy of Sciences carried out in collaboration with the US National Academy of Sciences with the financial support of external partner”, PAN.BFB.S.BWZ.329.022.2023. The work of PB was also supported by the Volkswagen Foundation under the special stipend No. 9D154. The work of CO was funded by the Science Committee of the Ministry of Science and Higher Education of the Republic of Kazakhstan, Grant No. BR20280974. This research is funded by the Science Committee of the Ministry of Science and Higher Education of the Republic of Kazakhstan, Grant No. BR21881880. This work has made use of data from the European Space Agency (ESA) mission Gaia (https://www.cosmos.esa.int/gaia), processed by the Gaia Data Processing and Analysis Consortium (DPAC, https://www.cosmos.esa.int/web/gaia/dpac/consortium). Funding for the DPAC has been provided by national institutions, in particular the institutions participating in the Gaia Multilateral Agreement. The use of TOPCAT, an interactive graphical viewer and editor for tabular data (Taylor 2005), is acknowledged.

References

Appendix A: Calculation of probabilities of membership to the group of Cr 135 and UBC 7

We adjusted our approach to selecting probable members of a system based on data obtained in our previous investigation. That is, in bluePaper I, we found the number of stars in a close vicinity to the considered pair of OCs with space parameters suitable to one cluster and kinematic parameters of another cluster. These stars were not included in the previous strict selection. However, we interpreted these stars as being probable members of a system and further computed their joint probability based on parallax and their joint kinematic probability by taking into account the parameters of both clusters.

A.1. Kinematic probability

As in bluePaper I, we found a halo of stars around the Cr 135 and UBC 7 pair with a location and kinematics combination that cannot be attributed to only one of the two clusters but rather to their system. This is why for the selection of space of coordinates and velocities, we did not distinguish between the two clusters but rather chose stars with both velocities and parallaxes sufficiently close to those of either of the two clusters: Cr 135 or UBC 7. We describe the details of the selection process below.

Provided that our intention was to search for coeval stars probably related to the pair of clusters in an extended region of sky, the kinematical selection principles needed to be respectively changed from 2D to 3D velocities to account for the projection effect. We used the modification of the convergent point method described by Röser et al. (2011). Of the initial dataset including 116975 stars, radial velocities were available only for 4575 sources. Thus, we selected probable members related to the system based on the components of their tangential velocity.

We expected that the members related to Cr 135 and UBC 7 follow certain distributions with respect to the mean 6D phase-space parameters cited in Table 1 (based on Gaia DR3). The respective space velocities were obtained as

$$\begin{array}{c}\hfill \left[\begin{array}{c}{U}^i\\ V^i\\ W^i\end{array}\right]=V_r^i·\left[\begin{array}{c}cos{l}^{i}cos{b}^i\\ sin{l}^{i}cos{b}^i\\ sin{b}^i\end{array}\right]+\frac{\kappa {\mu }_l^i}{{\varpi }^i}·\left[\begin{array}{c}-sin{l}^i\\ cos{l}^i\\ 0\end{array}\right]+\frac{\kappa {\mu }_b^i}{{\varpi }^i}·\left[\begin{array}{c}-cos{l}^{i}sin{b}^i\\ -sin{l}^{i}sin{b}^i\\ cos{b}^i\\ \end{array}\right],\end{array}$$(A.1)

where κ = 4.74047 is the transformation factor from 1 mas/yr at 1 kpc to 1 km/s; the indices are i = 1 for Cr 135 and i = 2 for UBC 7; the mean values for coordinates, parallaxes, proper motions, and radial velocities are from Table 1.

The kinematic probability $P_{\text{kin}}^k$ for a k^th star to belong to the investigated system is defined as

$$\begin{array}{c}\hfill P_{\text{kin}}^k=max\{\begin{array}{c}exp\{-\frac{1}{4}[{\left(\frac{{\mu }_l^k-{\mu }_l^{ex{p}_1^k}}{{\epsilon }_{{\mu }_l}}\right)}^2+{\left(\frac{{\mu }_b^k-{\mu }_b^{ex{p}_1^k}}{{\epsilon }_{{\mu }_b}}\right)}^2]\},\hfill \\ exp\{-\frac{1}{4}[{\left(\frac{{\mu }_l^k-{\mu }_l^{ex{p}_2^k}}{{\epsilon }_{{\mu }_l}}\right)}^2+{\left(\frac{{\mu }_b^k-{\mu }_b^{ex{p}_2^k}}{{\epsilon }_{{\mu }_b}}\right)}^2]\}\hfill \end{array}\end{array}$$(A.2)

where

$$\begin{array}{cc}\hfill {\mu }_l^{ex{p}_i^k}=& (-sin{l}_k·U_i+cos{l}_k·V_i)/(\kappa /{\varpi }_k);\hfill \\ \hfill {\mu }_b^{ex{p}_i^k}=& (-cos{l}_{k}sin{b}_k·U_i-sin{l}_{k}sin{b}_k·V_i+cos{b}_k·W_i)/\hfill \\ \hfill & (\kappa /{\varpi }_k)\hfill \end{array}$$(A.3)

are values of proper motion in (l, b) for given coordinates l_k, b_k, and ϖ_k expected for the stars sharing space motion (U, V, W)₁, (U, V, W)₂ of Cr 135 and UBC 7, respectively. Scattering parameters of the proper motion ε_μl, ε_μb were set to 1.0 mas/yr (about 2 km/s in tangential velocity).

A.2. Parallax probability

In the tangential plane, we limited our search of members by the coordinates selected in Sec. 3, but we did not discriminate between the stars by their position within this scope. The approach for selection by parallaxes is as follows: We used the probability distribution based on the mean parallaxes of Cr 135 and UBC 7 from Table 1 and the notion of mean dispersion of parallaxes of their core members. We intentionally restricted this value to avoid contamination with peripheral members of neighbouring clusters (primarily NGC 2451B). The resulting formula for the probability of membership based on parallax value and parallax error is as follows:

$$\begin{array}{c}\hfill P_{\varpi }^k=\{\begin{array}{c}1,3.23<{\varpi }_k<3.66\hfill \\ exp\{-\frac{1}{2}{\left(\frac{{\varpi }^k-3.23}{{\epsilon }_{\varpi }^k}\right)}^2\},{\varpi }_k\le 3.23\hfill \\ exp\{-\frac{1}{2}{\left(\frac{{\varpi }^k-3.66}{{\epsilon }_{\varpi }^k}\right)}^2\},{\varpi }_k\ge 3.66.\hfill \end{array}\end{array}$$(A.4)

The selected scheme assigns $P_{\varpi }^k=1$ to all stars in the line-of-sight between the centres of the pair of clusters and within about 7.5 pc (in naive approach d^k[pc] = 1000/ϖ_k) to and from each centre. The wings of the probability distribution were defined by

$$\begin{array}{c}\hfill {\epsilon }_{\varpi }^k=\sqrt{{\sigma }_{\varpi }^2+{\sigma }_{{\varpi }_k}^2},\end{array}$$

where ${\sigma }_{\varpi }^2=0.25$ mas is the approximate mean dispersion of parallaxes of the joint ensemble of reliable members of Cr 135 and UBC 7, and ${\sigma }_{{\varpi }_k}^2$ is the nominal Gaia EDR3 parallax error for each source.

A.3. Photometric probability

Isochrones for 50 Myr (according to bluePaper I) for Gaia DR3 passbands from Maíz Apellániz & Weiler (2018) were obtained from the Padova webserver CMD3.7 ³ based on the calculations by Bressan et al. (2012) for solar metallicity Z = 0.0152. A recent investigation by Spina et al. (2021) has confirmed close to solar abundance for Cr 135 ([M/H]= − 0.011). We took into account the relation between A_G/E(BP − RP) and (BP − RP)₀ of a star and A₀ (interstellar extinction at λ = 550 nm) from Gaia Collaboration (2018a) to adjust isochrones for photometric probability calculation at a fixed value of A₀ = 0.1^m, bluePaper I.

The photometric probability $P_{\text{ph}}^k$ was computed in the CMD (G, BP − RP) as a maximum of probabilities and was calculated separately for Cr 135 (cluster 1) and UBC 7 (cluster 2) as follows. The photometric probability was set as equal to one for the stars occupying the CMD in an area between the isochrones for single stars and for unresolved binaries of equal mass. For all the other stars, the photometric probability depended on a difference between magnitudes:

$$\begin{array}{c}\hfill \mathrm{\Delta }{G}^i=min(|G^i-G^b|,|G^i-G^f|),\end{array}$$(A.5)

where G^b, G^f are respectively brighter and fainter values of magnitude limiting the domain between the isochrones for single stars and for unresolved binaries of equal mass for a given colour BP − RP:

$$\begin{array}{c}\hfill P_{\text{ph}}^k=max(exp\{-\frac{1}{2}{\left(\frac{\mathrm{\Delta }{G}_1^k}{{\epsilon }_{G_1}^k}\right)}^2\},exp\{-\frac{1}{2}{\left(\frac{\mathrm{\Delta }{G}_2^k}{{\epsilon }_{G_2}^k}\right)}^2\}).\end{array}$$(A.6)

The value of ${\epsilon }_{\text{ph}}^i$ was defined from the mean dispersion of parallaxes of the joint ensemble of reliable members of Cr 135 and UBC 7, which is equal to 0.25 mas, and the individual photometric error ${\sigma }_G^i$ estimated by the flux and flux error in G.

All Tables

All Figures

	Fig. 1. Evolution of the distance between the centres of Cr 135 and UBC 7 for the DR2 (red) and DR3 (blue) initial conditions.
In the text

	Fig. 2. Comparison of observed positions of probable cluster members and of model stars of basic N-body simulations. The light green and light blue dots are the members of Cr 135 and UBC 7, respectively. The coordinate frame was obtained by parallel translation of the standard Galactic coordinate system.
In the text

	Fig. 3. Evolution of the tidal radius (upper panel) and bound particle number (bottom panel). The green line represents Cr 135, and the blue line indicates UBC 7.
In the text

	Fig. 4. Orbital evolution of Cr 135 and UBC 7. The colour-coding represents the time of integration. The filled and open black circles indicate the present positions of the clusters. Th black crosses indicate the positions where the clusters dissolve.
In the text

	Fig. 5. Hertzsprung-Russell diagrams of Cr 135 and UBC 7 (including all unbound stars). Red dots correspond to an age of 0 Myr, light blue to 50 Myr, green to 100 Myr, and black to 150 Myr. The crosses and the stars show two types of white dwarfs from the N-body simulation.
In the text

Fig. 6. A selection of Gaia DR3 data in the broad vicinity of the Cr 135 and UBC 7. Cr 135, UBC 7, and their closest neighbouring clusters according to Cantat-Gaudin et al. (2020), left panel. Refined sample of 60% of the members of Cr 135 and UBC 7 over the same sky area, right panel. The colour-coding indicates the MP according to Eq. (1). The circles in both panels represent the half-number radii of the clusters according to Cantat-Gaudin et al. (2020). See discussion in the text.

In the text

	Fig. 7. Components of the OC system. The probable members of Cr 135 are shown with red dots, probable members of UBC 7 are indicated with blue dots, and extended halo members are shown with yellow dots. The upper panel shows the lb Galactic coordinates sky plane; the bottom panel shows the ϖ,  b plane.
In the text

Fig. 8. Components of the OC system at the colour-magnitude diagram. Probable members of Cr 135 are shown with red dots, probable members of UBC 7 are indicated with blue dots (left panel), and probable members of the extended halo are shown with yellow dots (right panel). Grey dots are for background stars. The green line represents the isochrone for 50 Myr for single stars, and the magenta line is for the same isochrone for unresolved binary stars with identical components. The red line is for the zero-age main sequence.

In the text

Fig. 9. Fits of the King profile according to Eq. (5) for the observed radial density profiles, based on the confident members of the clusters. The curves present distributions including the residual background stars (black dots on the yellow lines), the distributions corrected for the background (circles with statistical error bars), and the fitted King profiles (solid red curves). The vertical dashed lines show the core and tidal radii (rc and rt). The colour bands show errors in the radii.

In the text

Fig. 10. Open clusters Cr 135 and UBC 7 in the tangential plane contacting the celestial sphere at the common centre of the mass of the clusters. The red dots show probable members of Cr 135, and blue dots show probable members of UBC 7. The circles show the core and tidal radii. Black dots show the centres of the clusters. The “+” sign shows the centre of mass, which is assumed to be at 292 pc from the Sun. Linear coordinates ξ and η are parallel to equatorial ones.

In the text

	Fig. 11. Cumulative number distribution of stars compared against the projected radius for Cr 135 (left panels) and for UBC 7 (right panels). The coloured lines represent different N-body randomisations. The upper panel presents the IMF sampling, and the bottom panels show the PV sampling. The shaded grey regions show the observational data with uncertainties.
In the text

	Fig. 12. Same as in Fig. 11 but compared against the residual magnitude of the proper motion.
In the text

Fig. 13. Outwards motion of the corona stars with respect to the pair of star clusters’ centre of mass. Smoothed magnitudes Vr of residual line-of-sight velocities of receding (approaching) stars from Gaia DR3 are shown by red (blue) lines. The mean residual tangential velocities in the picture plane Vt are shown by the green line. The grey line shows the mean tangential components of the runaway velocities obtained in N-body simulations for the Dehnen γ = 0 initial stellar profile with SFE = 0.05 after 50 Myr of evolution (Shukirgaliyev et al. 2021), where the shadow represents the uncertainty.

In the text