Removal of Ca²⁺, Mg²⁺, Ba²⁺ and Sr²⁺ from shale gas flowback water from the Sichuan Basin in China by chemical softening under the guidance of OLI Stream Analyzer: precipitation behaviors and optimization study

Technically, the global recoverable shale gas is 214.5 trillion cubic meters, meaning it could be an alternative to conventional natural gas in the future (UNCTAD 2018). The USA was the first country to develop shale gas, and in 2018, its total production of shale gas reached 625.8 billion cubic meters, accounting for 60% of the total domestic production of natural gas (EIA 2019). The US shale gas boom has fueled the ambitions of the oil and gas industry in China, with the first shale gas well fractured successfully in 2010, which thus marked the beginning of the shale gas industry in China. Currently, as one of four countries to commercialize shale gas, China is the second-largest shale gas producer in the world, producing 10.8 billion cubic meters in 2018.

In China, especially in the Longmaxi Formation in the Sichuan Basin, the water consumption of fracturing fluid used in a single well is about 16,000–40,000 m³. About 10–20% fracturing fluid is sent back to the ground in four weeks, while the remaining water will flow back as a gas-water mixture with daily water production of 1–10 m³/d at the later stage. Generally, the flowback water (FBW) contains a variable composition of pollutants, including high salinity, suspended solids, organic chemicals, heavy metals, and naturally occurring radioactive materials (NORMs), making it hard to purify from the perspective of wastewater treatment (Gregory et al. 2011).

Traditionally, particularly in the early stage of exploitation, recycling is the main disposal method (Butkovskyi et al. 2017). However, once fracturing halts, recycling is not feasible. In addition, environmental regulations might hinder the application of deep well injection (Estrada & Bhamidimarri 2016), and sometimes the transportation distance is long, resulting in a high cost of wastewater disposal. Hence, discharge has gained more attention as a feasible pathway for operators and academic researchers (Nasiri et al. 2017). In the process proposed for discharging, the key technology is assumed to be desalination, whose difficulties and cost increase with the elevation of total dissolved solids (TDS) (Onishi et al. 2017). It is reported that reverse osmosis (RO) works well when TDS content of the FBW is below 40,000 mg/L, and mechanical vapor recompression (MVR) performs better when TDS is over 40,000 mg/L (Gregory et al. 2011). Considering that the TDS in FBW from the Sichuan Basin is usually less than 40,000 mg/L, RO should be the core desalination process. However, the scaling risks cannot be avoided with RO. Potential scale, such as calcium, magnesium, barium, and strontium minerals, might hinder the achievement of maximum processing potential (Mohammadesmaeili et al. 2010). Therefore, it is crucial to set up the correct softening pretreatment in order to guarantee the efficiency and continuity of the desalination process.

Water softening, often called hardness removal, is usually performed by chemical precipitation, electrochemical precipitation, cation ion exchange, and nanofiltration (NF). In recent years, electrochemical precipitation has been studied a lot, including multistage electrochemical precipitation, paired electrolysis, and electrodialysis reversal (Lee et al. 2013; Yu et al. 2018; Sanjuán et al. 2019). However, their efficiency does not achieve a satisfactory level in actual wastewater, and energy consumption is also quite high. At the same time, Dong et al. (2016) applied cation ion exchange to remove over 98% Ca²⁺ in artificial water. Other methods, including modified ion exchange and combination processes like ultrasound-cation ion exchange, have also been employed to remove hardness (Entezari & Tahmasbi 2009; Mautner et al. 2019). However, the application of ion exchange in very hard water will increase backwash times, which may shrink the life cycle of resins. Modified NF has been studied to remove divalent ions Mg²⁺ (Liu et al. 2015; Das et al. 2019), which all neglect the membrane pollutants such as scale ions Ca²⁺, Ba²⁺, Sr²⁺, and organics in real complex water. Chemical precipitation is well-used because of its low cost and easy availability, especially in very hard and high TDS water. Hence, chemical precipitation was selected as the pretreatment for FBW disposal. Generally, hardness removal focuses on Ca²⁺ and Mg²⁺ and therefore neglects their group members Ba²⁺ and Sr²⁺, whose minerals also contribute to scaling (Esmaeilirad et al. 2015). Since the quality of FBW always varies with the flowback time, the concentration of hardness ions changes as well. Simulation software can simulate the change tendency of ions and provide guidance on economical reagent addition schemes. However, the complete simulation and process schemes of softening for variable water quality have rarely been reported.

Hence, the focus of this study was on the pretreatment optimization and precipitation behaviors of Ca²⁺, Mg²⁺, Ba²⁺ and Sr²⁺ based on typical FBW in the Sichuan Basin. OLI Stream Analyzer was used to simulate chemical softening and to provide a suggestion on pretreatment optimization. Jar tests were conducted to verify the accuracy of the simulation and to investigate the precipitation behaviors of Ca²⁺, Mg²⁺, Ba²⁺ and Sr²⁺. A chemical softening procedure was given.

The water samples were taken from a shale gas well (Sichuan Province, China) and collected in 25 L plastic containers until the containers were full. The containers were then sealed until the experiments were carried out to prevent CO₂ from escaping. Prior to each experiment, the plastic container was shaken for several minutes until the particles and sediments in it were re-suspended in the water. The key characteristics of this water are shown in Table 1. In order to eliminate the suspended solids and oil and grease, pre-coagulation was applied. Based on previous experiments, the process was as follows: coagulation with 50 mg/L ferric chloride (analytical grade, Cologne Chemicals co. Ltd., China) at a rate of 800 sec⁻¹ velocity gradient for 1 min; flocculation with 6 mg/L anionic polyacrylamide solution (10 million molecular weight, Cologne Chemicals Co. Ltd., China) at a rate of 40 sec⁻¹ velocity gradient for 20 min. After settling for 30 min, the supernatant was filtered through a 0.45 μm nylon filter membrane. The main cations and anions remained basically the same as in the raw water.

MINEQL + , PHREEQC and OLI Stream Analyzer are widely used for chemical equilibrium calculation. The difference between them is the activity-ionic strength models they apply. The activity-ionic strength models can be divided into ion-association models and ion-interaction models. The main ion-association models are listed in Table 2.

For all equations shown in Table 2, is ionic strength and a_i and b_i are parameters determined by the ion size, while A and B are temperature-related factors.

The most commonly used ion-interaction model is the semi-empirical Pitzer equation based on statistical mechanics, assuming all charged ions are independent individuals with electrostatic interaction and short-range interaction. Due to the complexity of the system, the calculation procedure is complicated to carry out manually. MINEQL+ uses the Davis equation to calculate the activity coefficient, which means the software is available when the ionic strength does not exceed 0.5 molality. Only below the limits can the model give accurate results, but clearly, the ionic strength of the water samples exceeds the limit. Thus, it was not suitable for this study. PHREEQC, a hydrogeochemistry calculation software developed by the United States Geological Survey (USGS), was embedded in the Pitzer model in version 3.0, making it possible to apply in a high ionic strength solution. The OLI Stream Analyzer software, developed by OLI Systems Inc., uses the Debye-Hückel model and the Pitzer model and is widely applied in the petrochemical industry. In the Sichuan Basin, the TDS of FBW is lower than 40,000 mg/L, corresponding to the ionic strength of less than 1 M, thus it is not possible to apply PHREEQC. Besides, TDS changes with flowback time, which requires a wider adaption of software. OLI Stream Analyzer software with both ion-association and ion-interaction model databases has already been proven to be a good fit for the oil and gas industry (Yang 2014). Hence, OLI Stream Analyzer software was selected to simulate the softening process in this study.

Batch studies were performed using a series of jar tests to simulate chemical softening treatment (i.e., rapid mix, sedimentation and filtration). An MY3000-6 series six-paddle standard jar-tester (Meiyu Instruments co., Ltd, Wuhan, China) was used to perform all jar-testing experiments. All experiments were conducted at room temperature (20 °C).

NaOH: precipitation behaviors of FBW were studied in 1 L beakers. The contents of the beaker were completely mixed with paddles. The jar test paddles were run at a rate of 800 sec⁻¹ velocity gradient to simulate a rapid mix. 1 mol/L NaOH (analytical grade, Cologne Chemicals Co., Ltd., China) was prepared. Based on the simulations, the dosages of NaOH from 200 to 800 mg/L at intervals of 100 mg/L were added and allowed to rapid mix for 90 min to ensure sufficient reaction. After rapid mixing, it took 30 min for the sediments to settle and the pH to stabilize. The supernatant was filtered through 0.45 μm nylon filter membranes and immediately acidified to pH < 3 with HNO₃ (1:1, analytical grade, Cologne Chemicals Co., Ltd., China). The water samples were immediately analyzed to avoid further reaction.

Na₂CO₃: 0.5 mol/L Na₂CO₃ (analytical grade, Cologne Chemicals Co., Ltd., China) was prepared. Based on the simulations, 500, 1,000, 2,000, 3,000, and 4,000 mg/L dosages of Na₂CO₃ were added. The procedure was the same as above.

Na₂SO₄: 1 mol/L Na₂SO₄ (analytical grade, Cologne Chemicals Co., Ltd., China) was prepared. Based on the simulations, 50, 100, 150, 200 and 300 mg/L dosages of Na₂SO₄ were added. The procedure was the same as above.

NaOH-Na₂CO₃ process: pretreatment optimization studies were also conducted in 1 L beakers. According to the comparison of simulation and experiments, 1 mol/L NaOH with 700 mg/L was the optimum addition. In view of this, 500, 1,000, 1,500, 2,000 and 3,000 mg/L dosage of Na₂CO₃ were then added. The procedure was the same as above.

Na₂SO₄-NaOH-Na₂CO₃ process: based on the simulations and experiments, 1 mol/L Na₂SO₄ with 100 mg/L and 1 mol/L NaOH with 700 mg/L were the optimum additions. In view of this, 500, 1,000, 1,500, 2,000 and 3,000 mg/L dosages of Na₂CO₃ were then applied. The procedure was the same as above.

The cations were analyzed using an inductively coupled plasma optical emission spectrometer (ICP-OES, ICAP 6300, Thermo). Precipitated solids were washed with deionized water (DI water), and dried in an electric oven for subsequent petrographic and morphological characterization. For X-ray diffraction (XRD) analysis, solids were affixed to sample slides and a PANalytical X'Pert Pro Diffractometer (Almelo, the Netherlands) using Cu Kα radiation was utilized to conduct the experiments. Scans were performed over a 2-theta range between 5° and 70°. Pattern analysis was performed using the computer software MDI Jade (Version 6.5). For the scanning electron microscopy (SEM) and energy-dispersive X-ray spectroscopy (EDS) analyses, samples were imaged with a ZEISS EV0 MA15 scanning electron microscope (Carl Zeiss Micrographics Co., Ltd). The typical accelerating voltage was 0.2–30 kV with continuous adjustment. The supporting EDS was used to conduct elemental analysis.

The solubility product constants (K_sp) of possible softening precipitates are listed in Table S1 in the Supplementary Material. Since the K_sp of magnesium hydroxide (Mg(OH)₂) is only 1.8 × 10⁻¹¹, Mg²⁺ would easily form Mg(OH)₂ under alkaline conditions. As the raw water was neutral to slightly acidic (pH 6.50), it was necessary to supply alkali to react with Mg²⁺. According to the simulation of the OLI Stream Analyzer (as shown in Figure 1(a)), dissolved Ca²⁺ and Sr²⁺ continued to decrease with NaOH dosage ranging from 0 to 300 mg/L. The more alkaline the solution was, the higher the amount of converted to (Silva et al. 2019), which would contribute to the formation of calcium carbonate (CaCO₃) and strontium carbonate (SrCO₃). Once the NaOH dosage increased to over 300 mg/L, the carbonate alkalinity in water had been largely consumed, and the main reaction shifted to the reaction of Mg²⁺ and OH⁻ corresponding to the increasing NaOH dosage. When the pH was raised to about 11 with the addition of 700 mg/L NaOH, dissolved Mg²⁺ decreased sharply to less than 0.1 mg/L. Also, a comparison of the XRD patterns of the dosages of 200 mg/L and 600 mg/L NaOH (Figure 2) verified that the low content of OH⁻ could not induce the precipitation of Mg²⁺, which was consistent with the result of EDS (Figure 3).

The simulation of the OLI Stream Analyzer showed that the concentration of Ba²⁺ decreased by about 28.56% initially but barely changed through the softening process. The precipitate formed was shown to be BaSO₄ rather than BaCO₃, which was independent of NaOH addition. That could be attributed to the supersaturation of BaSO₄ in the raw flowback water (saturation index, SI = 4). As a matter of fact, the solubility of BaSO₄ would decrease with the reduction of the temperature and pressure (Moghadasi et al. 2003), as the situation when FBW returns to the surface after the hydraulic fracturing, which could easily create the supersaturated state of BaSO₄ in the raw flowback water. The metastable state of supersaturated BaSO₄ in solution could even reach 6,700 times over the solubility limit (Monnin & Galinier 1988), but it may precipitate after a long period of storage and tests. This phenomenon was confirmed by the XRD and EDS analysis, in which the characteristic peaks of BaSO₄ were identified (Figure 2) and elemental Ba was detected (Figure 3), indicating the precipitation of BaSO₄.

At the dosage of 1,000 mg/L NaOH, the removal ratios of Ca²⁺, Mg²⁺, Ba²⁺, and Sr²⁺ were 59.83, 98.46, 40.59 and 32.79%, respectively. As shown in Figure 1(b), the removal ratios of Ca²⁺ and Sr²⁺ in the jar tests were 10.5% and 27.9% lower than the simulation results, while that of Ba²⁺ was 11.72% higher than the simulation, which might be attributed to the coprecipitation mechanism. In tests, because of the similar ionic radii, Ba²⁺ could substitute Ca²⁺ and Sr²⁺ in the lattice of CaCO₃ and SrCO₃ during precipitation. In Figure 2, the low peak of alstonite (CaBa(CO₃)₂) and paralstonite ((Ba, Sr)Ca(CO₃)₂) were identified, which indicated the coprecipitation of Ca²⁺, Ba²⁺ and Sr²⁺. The presence of soluble Mg²⁺ in the solution could inhibit the precipitation of CaCO₃ (Rioyo et al. 2018), which might be another explanation for the distinctions of the hardness ion concentration in the jar tests. Since perfect chemical equilibrium was hard to reach during the bench-scale tests, the removal ratios might be lower than expected. Except for CaBa(CO₃)₂ and (Ba, Sr)Ca(CO₃)₂ formed by coprecipitation mechanism in the jar tests, XRD and SEM/EDS analysis also proved that the simulations could adequately predict the main reaction in water softening (shown in Table 3). The reason for the distinction was that software for chemical equilibrium calculation such as OLI Stream Analyzer are all based on thermodynamic equilibrium, which ignores the kinetics of scale mineral precipitation such as coprecipitation (Bozau & van Berk 2013).

According to the simulation (Figure 4(a)), because of the low solubility of CaCO₃ and SrCO₃, Ca²⁺ and Sr²⁺ first decreased to nearly 13.96 mg/L and 3.51 mg/L with the addition of 1,000 mg/L Na₂CO₃. When the dosage changed from 1,000 to 3,000 mg/L, Ba²⁺ began to precipitate from the aqueous phase, leading to the decrease of Ba²⁺ from 100.93 mg/L to about 12.05 mg/L. Finally, Mg²⁺ slightly reduced by 15.99% even when 3,000 mg/L Na₂CO₃ was added, proving that Mg²⁺ could not be removed adequately by only supplying Na₂CO₃. Through the tests, the concentrations of Ca²⁺, Sr²⁺ and Ba²⁺ were reduced from 370.25 mg/L to 33.83 mg/L, 68.945 mg/L to 28.31 mg/L, and 141.9 mg/L to 8.47 mg/L, respectively. In contrast, only 34.48% of Mg²⁺ was removed. As shown in Figure 4(b), except Mg²⁺, tested Ca²⁺, Ba²⁺ and Sr²⁺ concentrations were higher than the values determined by the calculation. Magnesium carbonate minerals always exist as mixed phases, and Mg²⁺ could substitute for Ca²⁺, Sr²⁺ and Ba²⁺ into the lattice of corresponding carbonate (Thorstenson & Plummer 1977; Morse et al. 2007), which promoted the removal of Mg²⁺.

Another mechanism that could not be ignored was the salt effect. The high content of NaCl, more than 10,000 mg/L, would strengthen ion-ion interactions and promote solubility (Cao 2007). It might partly explain why Ca²⁺, Sr²⁺ and Ba²⁺ did not precipitate quantitatively, as demonstrated by the simulation. Coprecipitation would facilitate the precipitation of Mg²⁺ but delay that of other ions, while a salt effect mechanism restrained the removal of all the hardness ions, leading to the higher removal of Mg²⁺ and lower removal of Ca²⁺, Sr²⁺ and Ba²⁺.

The XRD patterns (Figure 5) and EDS analysis (Supplementary Material, Figure S1) showed that the Mg element appeared in the precipitates when 2,000 mg/L Na₂CO₃ was used, which was earlier than the dosage provided by the simulation (3,000 mg/L), since Mg²⁺ would coprecipitate with Ca²⁺ to form (Mg_0.03Ca_0.97)CO₃. Besides, the form of barium-containing precipitates was identified to be BaCa(CO₃)₂ instead of BaCO₃, verifying the coprecipitation mechanism again. Comparison of the solids formed between jar tests and simulations is presented in Table S3 in the Supplementary Material, which reveals that coprecipitation could be used to explain the inconsistency between the simulations and jar tests.

The significantly higher solubility (3.89 g/100 g, 20 °C) of barium hydroxide as compared to that of magnesium hydroxide (9.628 × 10⁻⁴ g/100 g, 20 °C) indicated that it was not a practical way to remove Ba²⁺ from raw water by adding alkali (Esmaeilirad et al. 2015). Since BaSO₄ is less soluble than other forms of barium, and its precipitation kinetics is quite fast (Zhang et al. 2014), it is feasible to apply Na₂SO₄ to precipitate Ba²⁺ in engineering.

According to the simulation results of the OLI Stream Analyzer (Figure 6(a)), Ba²⁺ rapidly decreased to about 0.5 mg/L at the dosage of 300 mg/L Na₂SO₄. Also, Sr²⁺ could combine with and to form solids out of the aqueous phase. However, BaSO₄ usually precipitates before the appearance of SrSO₄ in the solution when Ba²⁺ and Sr²⁺ coexist, due to the much slower precipitation kinetic of SrSO₄ when compared with BaSO₄ and Sr²⁺ being likely to coprecipitate with BaSO₄ (Zhang et al. 2014). Therefore, the main reaction after adding Na₂SO₄ was the BaSO₄ precipitation, leading to over 90% removal of Ba²⁺. Although the concentration of Ca²⁺ remained unchanged, it was slightly lower than that in raw water, which might be attributed to the precipitation of CaCO₃. The XRD patterns (Supplementary Material, Figure S2) verified the existence of CaCO₃, indicating the supersaturation state of CaCO₃ in raw water. Test results showed that a 20.05% decrease of Sr²⁺ still existed, which could be explained by the incorporation of Sr²⁺ into the lattice of BaSO₄. Sr²⁺ was prone to coprecipitate with BaSO₄, as evidenced by the characteristic peaks of Ba_0.75Sr_0.25SO₄ in the XRD pattern. Nevertheless, as can be seen in Figure 6(b), the removal ratio of Ba²⁺ could reach 99.28% when only 300 mg/L dosage of Na₂SO₄ was used, verifying the high efficiency of Na₂SO₄ to remove Ba²⁺ even with complicated water chemistry in real FBW. The removal ratios of Ca²⁺, Mg²⁺, Ba²⁺, and Sr²⁺ were 12.09, 6.69 99.28 and 20.05%, respectively.

Since both tests and simulation verified that the precipitation of Ba²⁺ would not begin until the residual concentrations of Ca²⁺ and Sr²⁺ were negligible, the dosage of Na₂CO₃ should be high enough to induce the precipitation of Ba²⁺. A dosage of up to 4,000 mg/L Na₂CO₃ yielded only 66.81% removal. By contrast, the removal ratio of Ba²⁺ under the application of 300 mg/L NaSO₄ could reach 99.28%, which validated the efficiency of adding NaSO₄ to remove Ba²⁺ in FBW.

Since the individual additions could not achieve satisfactory removal efficiency for all hardness ions, as mentioned above, a combinated process was considered. A NaOH-Na₂CO₃ process is often used to induce the precipitation of hardness ions in groundwater (Rioyo et al. 2018). Based on the results in Figure 1, 700 mg/L NaOH was selected as the dosage in the first stage because of the 87.23% removal ratio of Mg²⁺. This section focuses on the selection of Na₂CO₃ dosage for the second stage. As shown in Figure 7, the removal ratios of Ca²⁺, Mg²⁺ and Sr²⁺ reached about 97.27, 77 and 65.04%, respectively, when 3,000 mg/L Na₂CO₃ was added after 700 mg/L NaOH. The removal efficiency of Ca²⁺ and Sr²⁺ was 6% higher than with Na₂CO₃ alone. Previous study has shown that the NaOH-Na₂CO₃ process is an effective way to remove calcium and magnesium hardness (Wang et al. 2019). However, Ba²⁺ has always been ignored, which might be a cause for BaSO₄ scaling, especially in the oil and gas industry (Bageri et al. 2017). In fact, even when the Na₂CO₃ dosage was raised to 3,000 mg/L, residual Ba²⁺ was still about 42.5 mg/L. The low removal ratio of Ba²⁺ suggested that the NaOH-Na₂CO₃ process was not suitable for FBW with a high concentration of Ba²⁺.

As mentioned above, since 300 mg/L Na₂SO₄ could remove 99.28% of Ba²⁺ in 90 min in the individual process, the Na₂SO₄-NaOH-Na₂CO₃ process was applied. Li (2011) preliminarily verified the feasibility of dosing sulfate and carbonate in the chemical softening process, but there was no specific method in the study. Considering that Ba²⁺ could also be partly removed by coprecipitating with Ca²⁺ and Sr²⁺ through the addition of NaOH and Na₂CO₃, the dosage of 100 mg/L Na₂SO₄ was selected. When the dosage of 2,000 mg/L Na₂CO₃ was added, the concentrations of Ca²⁺, Mg²⁺, Ba²⁺, and Sr²⁺ were reduced to 14.7, 15, 4.43, and 26.1 mg/L, respectively (Figure 8). Through the process, the removal ratio of Ba²⁺ was 31.84% higher than that obtained by the NaOH-Na₂CO₃ process. Further, the market price of Na₂SO₄ is lower than that of Na₂CO₃, meaning that the treatment was both effective and economical.

The scaling tendency of the FBW was predicted by the Reverse Osmosis System Analysis (ROSA) software (Version 9.0). ROSA software was developed by the Dow Chemical Company to model and predict the performance of RO systems, which is based on empirical data generated by testing FILMTEC^TM products at different conditions (Boulahfa et al. 2019). As shown in Table 4, the scaling tendency of RO was predicted with FBW as inlet water at various recovery rates. In the table, the Langelier saturation index (LSI) and Stiff & David stability index (S&DSI) were used to evaluate the scaling possibility of CaCO₃. If they are greater than 0, the scaling risk of CaCO₃ exists. For other mineral scale, there is a risk of scaling if the saturation percentage exceeds 100%. With different recovery rates, LSI was greater than 0, which suggested the possibility of the precipitation of CaCO₃. However, more remarkable was the fact that even when BaSO₄ saturation far exceeded 100%, the recovery rate was only 10%.

In the operation of RO, pH regulation is often used to control the precipitation of CaCO₃, Mg(OH)₂, etc. However, it is not effective to prevent the scaling of BaSO₄, which can cause flux decline and membrane damage in RO (Boerlage et al. 1999). Thus, it is essential to remove Ba²⁺ in the softening process as a pretreatment to desalinate the FBW for the discharge process. After the softening treatment of the Na₂SO₄-NaOH-Na₂CO₃ process with the dosage of 100 mg/L (1 mol/L), 700 mg/L (1 mol/L), and 2,000 mg/L (0.5 mol/L), no scaling tendency was predicted.

This study aimed to provide a suitable chemical softening method to prevent mineral scaling in RO desalination during shale gas flowback water treatment prior to discharge. The results showed that high concentrations of Ba²⁺ in flowback water need to be removed during the softening process, since they cause significant BaSO₄ scaling. The conventional NaOH-Na₂CO₃ process could precipitate Ca²⁺, Mg²⁺ and Sr²⁺ effectively, but the removal ratio of Ba²⁺ was limited unless a high dosage of Na₂CO₃ was applied. A small dosage of Na₂SO₄ could lead to a remarkable reduction of Ba²⁺ and also save on Na₂CO₃. OLI Stream Analyzer could predict the precipitation behaviors of Ca²⁺, Mg²⁺, Ba² and Sr²⁺ in the chemical softening process. This study presents an effective chemical softening method and might provide guidance for shale gas flowback water treatment, which would contribute to suitable shale flowback water solutions and eventual sustainable shale gas development in China.

This work was financially supported by a grant (2015SZ0007) from the Science & Technology Department of Sichuan Province, China, and the Science and Technology Cooperation Project of the CNPC-SWPU Innovation Alliance. The authors gratefully acknowledge PetroChina Company Limited and its branch company PetroChina Southwest Oil & Gasfield Company for financial support and assistance in the fieldwork and simulation.

All relevant data are included in the paper and its Supplementary Information.