Yanping Duan ,Pengfei Wang ,Wenge Yang ,Xia Zhao ,Hong Hao, *,Ruijie Wu ,Jie Huang

1 School of Chemical Engineering,Northwest University,Xi’an 710069,China

2 Luoyang JALON Micro-Nano New Materials Co.,Ltd,Luoyang 471000,China

Keywords:Kinetic hydrate inhibitors Synthesis Poly(N-vinylcaprolactam-co-butyl methacrylate)Natural gas Hydrate Computer simulation

ABSTRACT Natural gas hydrate inhibitor has been serving the oil and gas industry for many years.The development and search for new inhibitors remain the focus of research.In this study,the solution polymerization method was employed to prepare poly(N-vinyl caprolactam-co-butyl methacrylate) (P(VCap-BMA)),as a new kinetic hydrate inhibitor (KHI).The inhibition properties of P(VCap-BMA) were investigated by tetrahydrofuran(THF)hydrate testing and natural gas hydrate forming and compared with the commercial KHIs.The experiment showed that PVCap performed better than copolymer P(VCap-BMA).However,low doses of methanol or ethylene glycol are compounded with KHIs.The compounding inhibitors show a synergistic inhibitory effect.More interesting is the P(VCap-BMA)-methanol system has a better inhibitory effect than the PVCap-methanol system.1% P(VCap-BMA)+5% methanol presented the best inhibiting performance at subcooling 10.3 °C,the induction time of natural gas hydrate was 445 min.Finally,the interaction between water and several dimeric inhibitors compared by natural bond orbital(NBO)analyses and density functional theory(DFT)indicated that inhibitor molecules were able to form the hydrogen bond with the water molecules,which result in gas hydrate inhibition.These exciting properties make the P(VCap-BMA) compound hydrate inhibitor promising candidates for numerous applications in the petrochemical industry.

1.Introduction

Gas hydrates are crystalline compounds composed of water cages and gas molecules,which are stabilized at a low temperature and high pressure,and often encountered in deepwater offshore operations [1-3].Gas hydrates present a very real challenge to the gas and oil industry.For plugging of the oil and natural gas transmission pipelines by the hydrates can cause financial loss and safety risks [4].Traditionally,a large quantity (up to 50%(mass))of alcohol and glycol was injected into production streams to prevent hydrate formation [5-7].The cost of chemicals and clean-up has urged an interest in what is called low dosage hydrate inhibitors(LDHIs).Kinetic hydrate inhibitors(KHIs)are one type of LDHIs,can inhibit the formation of natural gas hydrate under a low concentration of 0.5%-3%.KHIs are generally water-soluble polymers,for example,poly(vinyl caprolactam) (PVCap),poly(vinyl pyrrolidone) (PVP),antifreeze protein-based polymers,polyaspartamides,and their copolymers [8,9],even natural materials,such as polysaccharides and whey protein [10,11].Since most KHIs polymers have both hydrophobic and hydrophilic groups,they could regard as surfactants[12,13].More importantly,gas hydrate inhibitors can form a comparatively stronger hydrogen bond with water molecules of the clathrate hydrate to break the hydrogenbonded network of the clathrate structure [14].

Although KHIs would reduce the growth rates of hydrate crystal by binding to growth sites,the hydrate fractional number will ultimately become higher in which the hydrate particles would deposit and clog the hydrocarbon flow [15].Some studies have been made to combine kinetic inhibitors with thermodynamic inhibitors to achieve synergistic effects [16,17].For example,synergism between ethylene glycol(MEG)and PVP has demonstrated by mean of the increase in the induction time and the improvement of sub-cooling effects when 4%of MEG was added into a mixture solution of 2% (mass) PVCap and PVP (an mole ratio of 1:1)[18].Renatoet al.research showed that ethanol has good synergy with KHIs.When 10%of ethanol is added,the induction time of the combination of the two tested inhibitors increases significantly[16].

Several mechanisms have been put forward to explain the inhibition of gas hydrates using kinetic hydrate inhibitors.KHIs are used to modify the local structure of the water molecule and cause the inhibition of hydrate nucleation.In the presence of KHIs,prevent the local organization of guest and host molecules to enhance the nucleation barriers.KHIs can attach to the nucleated surface and stop further hydrate growth [19,20].Moreover,recently it is reported another effect of KHI in the literature as it can cause hydrate dissociation inside the hydrate phase boundary [21].At the same time,studying the mechanism of action of inhibitors through computational simulation is also the focus of attention.For example,hydrate formation or decomposition kinetic has researched in the presence of different KHIs,and the property of the interaction between an inhibitor and liquid water has been studied as a two-stage computational from initially density functional theory (DFT)calculations to molecular dynamic(MD)simulation [22].Kuznetsova [23] has revealed that placing PVP into a methane-water system and rapidly diffuses to the methanewater interface and stays there,moving alongside the surface throughout the imitation.Johnet al.[24] has discovered that PVCap containing two monomer repeat units give a structural insight into the mode of interaction among low dosage hydrate inhibitors and water,shown that this simple compound represents an appropriate model for commercial KHIs.Furthermore,the quantum chemical calculation could play a significant role in understanding hydrogen bond interaction [14].

In this paper,a new kinetic inhibitor P(VCap-BMA)was synthesized and compared with commercial KHI such as PVCap and Inhibex 501.According to experimental results,P(VCap-BMA) has weaker inhibition properties than PVCap in the same condition.While added low dosage ethylene glycol or methanol into KHI,the complex hydrate inhibitor exhibits synergism,and improve the inhibition performance of P(VCap-BMA).Lastly,using NBO analyses and density functional theory (DFT) calculated the interaction energy between water molecules and inhibitor molecules to further understand the mechanism of action of inhibitors.

2.Experimental

2.1.Materials

Azobisisobutyronitrile (AIBN,99%) was supplied by Tianjin Guangfu Fine Research Institute.N-Vinylcaprolactam (Vcap,98%)was purchased from Sigma-Aldrich.Tetrahydrofuran (THF,≥99.0%) was purchased from Tianjin Fuyu Fine Chemical Co.,Ltd.N-hexane (95.0%) was purchased from Tianjin jizhun Chemical Co.,Ltd.Anhydrous diethyl ether (≥99.5%) was purchased from Sichuan Xilong Chemical Industry Co.,Ltd.N-butyl methacrylate(BMA,99.5%) was purchased from Tianjin Kemiou Chemical Reagent Co.,Ltd.The natural gas (CH492.05 mol%) used was purchased from Beijing Tian youshun Gas Co.,Ltd,the composition shown in Table 1.Commercial inhibitor inhibex 501 (2-ethylene glycol monobutyl ether solution of poly(N-vinyl acetamide-N-vinyl caprolactam)) was purchased from the company of Ashland Specialty Chemicals.PVCap was synthesized by our laboratory(Mn=8792).

2.2.Fabrication of the P(Vcap-BMA) copolymer

The copolymer P(Vcap-BMA) was prepared by the solution polymerization method.22.00 g of Vcap,2.00 g of BMA,0.05 g of AIBN and 50 ml distilled water were added into a 100 ml reactionflask.The reaction mixture was purged with nitrogen for 30 min and then heated to 75 °C under stirring at 400 r·min-1,reacted for 6 h.Pour the solution after the reaction into a 100 ml beaker and standing at 45 °C for 12 h,removed the supernatant liquid,and obtained a crude product.Finally,the crude product was dissolved in THF,extracted withn-hexane,and washed with anhydrous diethyl ether.After two precipitation cycles and removing the solvent by a rotary evaporator,the copolymers were dried out a constant weight in a vacuum oven at 37 °C under 0.03 MPa.Structures of PVCap and P(VCap-BMA) shown in Fig.1.

Table 1 Composition of synthetic natural gas used in tests

2.3.Characterization

Fourier Transform Infrared (FT-IR) spectra were recorded from 400 cm-1to 4000 cm-1on a Nicolet MX-1 FT-IR (Nicolet,USA)the spectrometer,using the KBr pellet technique.1H nuclear magnetic resonance (1H NMR) spectrograms of polymers were determined by a Brucker AVANCE 400 NMR instrument (Brucker,Germany) at room temperature,deuterated chloroform (CDCl3)used as the solvent and tetramethylsilane (TMS) as the internal standard.The molecular weight of polymers,using polystyrene samples as molecular-weight standards,was determined by GPC-1500gel permeation chromatography (GPC) (Waters,USA) a 28°C.Thermogravimetric (TG) performance was examined by a Netzsch TG 209 (NETZSCH,Germany) machine in the nitrogen atmosphere at a scanning rate of 10 °C·min-1from room temperature to 900 °C.

The cloud point (TCl) measurement was as follows:0.5 g of P(VCap-BMA) dissolved in 100 ml of deionized water.5 ml of the solution was placed in a transparent tube and heated it slowly until there was a haze in the container.TheTCldecided at the point when the first sign of fog in the solution was observed.Every test repeated three.

The results of the1H NMR analysis,IR measurement,and GPC for P(VCap-BMA) are listed in Table 2.

TG analysis of P(VCap-BMA) copolymer as shown in Fig.2.At 480 °C,96% of mass loss observed to occur owing to the thermal decomposition of the copolymer chain.The comparatively higher thermal stability of the P(VCap-BMA)could find at the temperature below 300 °C.

Fig.1. Structures of PVCap (a) and P(VCap-BMA) (b).

Table 2 The characteristics of P(VCap-BMA) used in this work

Fig.2. TG and DTG curve of P(VCap-BMA) copolymer.

2.4.THF hydrate inhibition testing

A water solution of THF was used as a model for natural gas because of the similar thermal and mechanical properties between the THF hydrate (structure II) and methane hydrate (structure I).THF hydrates crystals formed from 18.9% (approx.molar ratio 1:17) THF aqueous solution at atmospheric pressure,and its equilibrium temperature is 4.4 °C [25].

The experimental procedure was as follows:(1)A total of 80 ml water solutions of the THF was added into a 100 ml glass flask.(2)The experimental chemical was dissolved in this solution to give the desired concentration;for example,0.8 g of P(VCap-BMA) in 80 ml of this solution gave a 1% (mass) solution of P(VCap-BMA).(3)The round flask was in a cooling bath preset to set the temperature at(0±0.1)°C,which represents about 4.4°C subcooling.The degree of subcooling is the difference between the hydrate equilibrium temperature and the current experimental temperature under the test pressure.(4) The solution was slightly stirred manually with a glass rod without touching the glass flask walls.The formation time of THF hydrates was recorded to assess the KHIs inhibition performance.five parallel experiments were performed for each test solution,and the average induction time was used as the ultimate induction time of the inhibitor solution at this concentration.Ethylene glycol or methanol with the KHI was added into the THF solution when the synergistic effect of the complex hydrate inhibitor was investigated.

2.5.High-pressure kinetics inhibition testing.

A constant temperature method can measure natural gas hydrate induction time,near actual conditions.Fig.3 illustrates the experimental apparatus used in the assessment of the inhibition behavior of the KHIs.The details of this equipment and testing process were described in previous studies [12,26],with the central part being a cylindrical high pressure(up to 50 MPa)transparent sapphire cell with a sufficient volume of 40 cm3.The accuracy of pressure prob is±0.02 MPa,the accuracy of temperature prob is±0.1°C.The 10 ml distilled water and the selected inhibitor placed in the sapphire cell,the temperature of the air bath was set to 2°C,the pressure of the sapphire reactor was about 5 MPa.

2.6.Computational simulations

Understanding the interaction between inhibitors and water is very important to explore the mechanism of action of inhibitors.In this part,we take P(VCap),polyvinylpiperidone (PVPip),PVP and P(VCap-BMA) as the research objects.Gauss View 5.08 software was used to establish the interaction model between their dimer molecules and water molecules.The 6-311G (d,p) basis set under the DFT/M06-2X method in Gaussian 09 [27] was used to geometrically optimize the model.Charge transfer and interaction energy between water molecules and dimer inhibitor molecules were analyzed by natural bond orbital (NBO) analysis and density functional theory (DFT).

Fig.4 is the optimized energy configuration structure diagram of the dimer inhibitor molecule and water molecule model,in which P(VCap-BMA) a molecule and a water molecule established two models,P(VCap-BMA)1-H2O refers to the interaction between the oxygen atom on the amide ring of P(VCap-BMA)and the hydrogen atom in water,named as the model I.P(VCap-BMA)2-H2O is the interaction of the ester carbonyl oxygen atom of P(VCap-BMA)with a hydrogen atom in water,named as model II.The interaction model of PVCap-H2O,PVPip-H2O,PVP-H2O and H2O-H2O is named model III,IV,V,VI respectively.

To correctly characterize the interactions between the inhibitor molecule and a water molecule,we analyzed the calculation of the change of charge distribution in the stable configuration with the NBO method,a standard option in Gaussian 09 [28].The charge distribution on the atomic sites deduced using the NBO method at the DFT/M06-2X level[29].The interaction energy(ΔE)between the inhibitor and the water mixing system could be calculated using Eq.(1).The interaction energies were calculated at the level M06-2X/6-311 g(d,p)basis of theory and are corrected for basis set superposition error (BSSE) with the counterpoise method [30].

in whichEijis counterpoise corrected energy of the system withiandjmolecules inside,andEiandEjare the optimized energy ofimolecule andjmolecule,respectively.

3.Results and Discussion

3.1.Effect of polymer structure on the THF hydrate

Fig.3. A schematic diagram of an experimental device.

To gain insight into the influence of the mass ratio of VCap and BMA on the inhibition performance,polymerizations are executed with different mass ratios of VCap to BMA.The induction time of THF hydrate of P(VCap-BMA) with different monomer mass ratios at 4.4 °C subcooling is listed in Table 3.

P(VCap-BMA)was prepared by radical solution polymerization.The mass ratio of VCap and BMA was 11:1,the dosage of the initiator was 2.5%(mass)of the total mass,and the temperature was 75°C.The measured THF hydrates induction time was 138 min.As this mass ratio is increased from 3:1 to 11:1,induction time increases,while 13:1,induction time decreases slightly.Hydrophobic side-chain with ester groups in P(VCap-BMA)formed a stronger steric effect between liquid water and hydrate nuclei,which makes hydrate nuclei growth more difficult.However,having more hydrophobic groups led to the lower solubility of P(VCap-BMA)in water and the inhibition to be less effective.

At the formation of THF hydrates,when hydrates were visually observed,that time expresses the hydrate induction time [31].Table 4 showed the induction time of the THF hydrates with different concentrations of P(VCap-BMA) or PVCap,respectively,and their cloud points.The induction time of THF hydrate increased with the increase of the inhibitors content in THF solution.This shows that the interaction between PVCap and P(VCap-BMA) and hydrates is enhanced,resulting in a further reduction of hydrogen-bond interaction between water and the hydrate surface,which were more unfavorable for hydrate formation.P(VCap-BMA)had a worse THF hydrate inhibition performance than PVCap at the same concentration of weight.In contrast to PVCap(Tc1=38),the presence of P(VCap-BMA) has a lower cloud point of 36 °C,due to the presence of the acrylate segment contained in P(VCap-BMA) enhanced the overall hydrophobicity.Further mechanism analysis will be discussed in the calculation simulation in Section 3.4.

3.2.Effect of KHIs and THIs compound inhibitor on THF hydrate formation

High dosages of methanol or glycol are used as thermodynamic inhibitors (THIs) to prevent hydrate formation and pipe blockage.However,the cost of operation and recovery is relatively high.The combination of kinetic inhibitors and low-dose thermodynamic inhibitors can achieve the purpose of synergistic increase,while reducing the amount of thermodynamic inhibitors [32].In a gas field,if the 30% (mass) of methanol is required to avoid gas hydrate problems,low doses of KHI could potentially replace 20%methanol [33].

In this work,we choose methanol and ethylene glycol for compounding with KHIs,because they are widely used thermodynamic hydrate inhibitors (THIs) in the gas fields.Table 5 shows the hydrate induction time of THF solutions containing different concentrations of KHIs and methanol or ethylene glycol when a subcooling 4.4 °C or 10.4 °C.

As shown in Table 5,compared to the effect of individual inhibitor (Table 4),behaviors of hydrate formation obviously delayed in the KHIs+THIs system when the KHI was added to 1%-5%methanol or ethylene glycol solution.At the same condition,KHI+methanol system was better THF inhibition performance than the KHI+ethylene glycol system.For 1% PVCap and 1% P(VCap-BMA),When the mass amount of methanol added is greater than or equal to 3%,P(VCap-BMA)+methanol system had a better THF inhibition performance than PVCap+methanol system.Compounding of 1% P(VCap-BMA) and 5% methanol shows the best inhibition performance.The THF hydrates formation is delayed for over 3300 minutes.For 3% P(VCap-BMA),added 3% or 5% THIs.Delayed behaviors of hydrate formation were more obvious.The 3% P(VCap-BMA)+5% THI system has a good THF hydrate inhibition performance even at a subcooling degree of 10.4 °C.Which probably because the solubility of P(VCap-BMA) was improved by the addition of methanol or ethylene glycol,at the same time the hydrophobic polymer side-chains might reduce the affinity of the THF (guest) in the water phase toward the hydrate crystals,lead to even more efficient inhibition of hydrate.

3.3.Effect of different inhibitors system on gas hydrate formation

The natural gas hydrate induction time of several kinetic inhibitors close to the actual condition was measured by a visible transparent sapphire tester.The recorded pressure-time relationship curves for different inhibitor systems are shown in Fig.5.The inhibition test results on the gas are listed in Table 6,including the corresponding subcooling.

Fig.4. The optimized energy configuration structures of inhibitor-water (I-V) and water-water (VI) complexes.The dotted line is a hydrogen bond.

Table 3 The effect of the mass ratio of monomer on THF hydrates average induction time at subcooling of approximately 4.4 °C.

Table 4 THF hydrates induction time for two inhibitors and TCl values at 0.5% (mass)

The typical pressure variation usually undergoes gas dissolving(the first rapid drop stage),induction of nucleation,and the growth of hydrate (the second rapid drop stage) three stages for an enclosed hydrates formation system.The induction time of gas hydrate formation is the time between the beginning of the dissolution of gas in the water phase and the start of the second rapid descent phase [26].

As showen in Table 6,to the condition of a subcooling of about 10.4°C and under about 5000 kPa pressure.1% P(VCap-BMA)+5%methanol showed the best inhibiting effects and the induction time of natural gas hydrates is 445 min.Commercial inhibitor inhibex 501 with a mass fraction of 0.5% has a moderate inhibitory effect.For 1%P(VCap-BMA)+5%ethylene glycol system,the induction time was 114 min.While 1% P(VCap-BMA) and 1% PVCap showed a weak inhibiting effect and an induction time of 28 min and 81 min,respectively.For these test inhibitors,the inhibitionperformance of natural gas hydrates and THF hydrates showed the same trend.

Table 5 THF hydrate average induction time of different KHI+THI systems at subcooling 4.4 °C or 10.4 °C under atmospheric pressure

Fig.5. The time -pressure graphs acquired under the constant temperature procedure at 2 °C.

3.4.NBO and interaction energy analysis

Charge transfer plays a significant role in hydrogen bond formation.It stands for electron delocalization interaction from the occupied molecular orbital of one molecule to the unoccupied molecular orbital of another molecule [34].The charge transfers that calculated by using NBO for P(VCap-BMA)1-H2O complex (I),P(VCap-BMA)2-H2O complex of (II),PVCap-H2O complex (III),PVPip-H2O complex (IV),PVP-H2O complex (V),and water-water(VI) are presented in Table 7.The charge distribution of atoms was listed in Table S1 and Table S2 in Supporting Material.The hydrogen bond of the main atoms is marked out by the dotted line in Fig.4.

Table 6 Induction time of gas hydrate formation in the presence of different kinds of inhibitors at subcooling of approximately 10.4 °C

Table 7 Summary of natural charge analysis between the inhibitor and water complex

From Table 7,for the model I,the net charge of a water molecule showed a decrease of 0.0113e,and a corresponding net charge of the inhibitor P(VCap-BMA) showed an increase by 0.0113e.For model II,by contrast,the net charge of water molecules decreased by 0.0065e,and the net charge of corresponding inhibitor P(VCap-BMA) increased by 0.0065e.These indicate that the interaction between the oxygen atom on the amide ring of P(VCap-BMA)and the hydrogen atom in the water was dominant.For III,IV,V,and VI,net charges of a water molecule decreased,and the net charges of the corresponding molecule of PVCap,PVPip,PVP,and H2O in turn increased.In a word,for models I,II,III,IV,V,and VI,the reduction in the charge of water molecules was the same as the increase in the complex inhibitor molecule.

As shown in Fig.4,from the optimized stable configuration,it is found that the distance between O(20) and H(51) atoms in the model I is 0.1857 nm,and the distance between O(50) and H(22)is 0.2648 nm.In model II,the spacing between O(14) and H(51)is 0.1889 nm,the spacing between O(50) and H(39) is 0.2340 nm.In model III,IV,V,and model VI,the shortest hydrogen bond distance (O···H) is (O20···H51) 0.1834 nm,(O7···H45)0.1889 nm,(O9···H3) 0.1840 nm,and (O1···H6) 0.1908 nm,respectively.The distance between atoms affects the amount of charge transferred between them to a certain extent.The longer the bond length,the lower the corresponding bond energy.From the two different bond lengths in the same configuration,it canbe seen that the intermolecular hydrogen bond between water and the inhibitor mainly occurs between the hydrogen of the water molecule and the oxygen atom of the inhibitor molecule.

Table 8 Interaction energy between one water molecule and 2-monomer inhibitor calculated at the M06-2X/6-311 g(d,p) level of density functional theory

The interactive energy between the water molecule and dimer inhibitor molecule is calculated using M062X/6-311 g (d,p) to study the relationship attractive or repulsive interaction between an inhibitor molecule and a water molecule.The negative interaction energy between numerators represents an attractive interaction,while the positive energy expresses a repulsive interaction.The calculated interaction energy is listed in Table 8.

For inhibitors studying in Table 8,the interaction energy of all inhibitor-water complex had negative values.The more negative values,the stronger the attractive interaction.For all models,the interaction energy from small to large is II ＜VI ＜IV ＜I ＜V ＜III.Except for model II,they are all more negative than the interaction energy in water molecules (-22.8419 J·mol-1).Consequently,all inhibitor molecules can form hydrogen bonds with the water molecule.The stronger a binding of the inhibitor is to water,and the more disruptive an inhibitor is to the structure of forming a hydrate nuclei [35-37],which causes gas hydrate inhibition.The interaction energy of the PVCap-water complex is-52.5100 kJ·mo l-1,which was negatively highest among the inhibitors.The interaction energy in model II and I are -13.3900 kJ·mol-1and -36.7 570 kJ·mol-1,respectively,smaller than -52.5100 kJ·mol-1in model III,indicating that the hydrophilicity of PVCap was better than P(VCap-BMA).For P(VCap-BMA),the smaller interaction energy to water indicated an impressive agreement between the calculated and experimental results.

4.Conclusions

In this paper,inhibitor P(VCap-BMA) is synthesized by copolymerization ofN-vinylcaprolactam and butyl methacrylate.Kinetic inhibition properties of P(VCap-BMA) and PVCap were evaluated by THF hydrate method and high-pressure sapphire test method.P(VCap-BMA) has a worse inhibition performance than PVCap as the hydrophobic groups of the P(VCap-BMA) side-chain led to a lower solubility in water and weaker hydrogen bonds with the hydrate surface.

However,when P(VCap-BMA) and PVCap are combined with low-dose methanol or ethylene glycol to make a compound inhibitor,the KHI+methanol system has a better synergistic effect than KHI+ethylene glycol system.And P(VCap-BMA)+methanol composite inhibitor has better inhibition performance than PVCap+methanol system.1% P(VCap-BMA)+5% methanol presented the best inhibiting results on the THF solutions,and the induction time of hydrates was over 3300 min.1% PVCap+5%methanol exhibited the induction time of hydrates was 887 min.Because the solubility of P(VCap-BMA)was improved with the help of alcohols,and the hydrophobic polymer side-chains and a stronger steric hindrance might reduce the affinity of the guest in the water phase for the hydrate crystals,lead to more efficient inhibition of hydrates.

The interaction relationship between the dimer inhibitor molecule and the water molecule is simply established through NBO and DFT.The results show that there are a charge transfer and mutual attraction between the water molecule and the inhibitor molecule.The interaction energy between PVCap molecules and water is the strongest.After the introduction of butyl methacrylate copolymerization,the ester group on the hydrophobic side chain repels the water molecule,the interaction between the oxygen atom on the amide ring and the water molecule is weakened,and the hydrophobicity of P(Vcap-BMA) is enhanced.

Declaration of Competing Interest

The authors declare that they have no known competing financial interests or personal relationships that could have appeared to influence the work reported in this paper.

Acknowledgements

We are grateful to Prof.Chang-Yu Sun(China University of Petroleum)for enlightening discussions and visual observation experiments.This work is supported by the Key Science and Technology Program of Shaanxi Province (2014K10-03).

Supplementary Material

Supplementary data to this article can be found online at https://doi.org/10.1016/j.cjche.2020.10.003.