Anatomy of the massive star-forming region S106. The [O I] 63 μm line observed with GREAT/SOFIA as a versatile diagnostic tool for the evolution of massive stars

The central area (40" x 40") of the bipolar nebula S106 was mapped in the [O I] line at 63.2 mu m (4.74 THz) with high angular (6") and spectral (0.24 MHz) resolution, using the GREAT heterodyne receiver on board SOFIA. The spatial and spectral emission distribution of [01] is compared to emission in the CO 16 -> 15, [C II] 158 mu m, and CO 11 -> 10 lines, mm-molecular lines, and continuum. The [01] emission is composed of several velocity components in the range from -30 to 25 km s(-1). The high-velocity blue and red-shifted emission (v = 30 to 9 km s(-1) and 8 to 25 km s(-1)) can be explained as arising from accelerated photodissociated gas associated with a dark lane close to the massive binary system S106 IR, and from shocks caused by the stellar wind and/or a disk envelope interaction. At velocities from 9 to 4 km s(-1) and from 0.5 to 8 km s(-1) line wings are observed in most of the lines that we attribute to cooling in photodissociation regions (PDRs) created by the ionizing radiation impinging on the cavity walls. The velocity range from -4 to 0.5 km s(-1) is dominated by emission from the clumpy molecular cloud, and the [O I], [C II], and high-J CO lines are excited in PDRs on clump surfaces that are illuminated by the central stars. Modelling the line emission in the different velocity ranges with the KOSMA-tau code constrains a radiation field chi of a few times 10(4) and densities n of a few times 10(4) cm(-3). Considering self-absorption of the [01] line results in higher densities (up to 106 cm 3) only for the gas component seen at high blue and red velocities. We thus confirm the scenario found in other studies that the emission of these lines can be explained by a two-phase PDR, but attribute the high-density gas to the high -velocity component only. The dark lane has a mass of similar to 275 M-circle dot and shows a velocity difference of similar to 1.4 km s(-1) along its projected length of similar to 1 pc, determined from (HCO+)-C-13 1 -> 0 mapping. Its nature depends on the geometry and can be interpreted as a massive accretion flow (infall rate of similar to 2.5 x 10(-4) M-circle dot yr(-1)), or the remains of it, linked to S106 IR/FIR. The most likely explanation is that the binary system is at a stage of its evolution where gas accretion is counteracted by the stellar winds and radiation, leading to the very complex observed spatial and kinematic emission distribution of the various tracers.