High-Frequency Pulsatile Parameterization Study for the Titania Ceramic Membrane Fouling Mitigation in Oily Wastewater Systems Using the Box–Behnken Response Surface Methodology

In this comprehensive study, a seven-channel ultrafiltration (UF) titania membrane was used to investigate the impact of the pulsatile cleaning process on the crossflow filtration system. Seventeen experimental runs were performed for different operating conditions with a transmembrane pressure (TMP) varying from 0.5 to 1.5 bar, a crossflow velocity (CFV) ranging from 0.5 to 1 m/s, and pulsatile parameters within an interval varying from 60 to 120 s with a duration of 0.8 s, and collecting membrane permeate flux and volume data. The optimized operating conditions revealed that a TMP of 1.5 bar, a CFV of 0.71 m/s, and a pulsatile cycle of 85 s were the best operating conditions to reach the highest steady permeability flux and volume of 302 LMH and 8.11 L, respectively. The UF ceramic membrane under the optimized inputs allowed for an oil-rejection ability of 99%. The Box–Behnken design (BBD) model was used to analyze the effect of crossflow operating conditions on the permeate flux and volume. The analysis of variance (ANOVA) indicated that the quadratic regression models were highly significant. At a 95% confidence interval, the optimum TMP significantly enhanced the flux and permeate volume simultaneously. The results also demonstrated a positive interaction between the TMP and the pulsatile process, enhancing the permeate flux with a slight impact on the permeate volume. At the same time, the interaction between the CFV and pulsatile flow improved the permeability and increased the permeate volume.


Introduction
Produced water (PW) is a significant byproduct of the oil and gas industry. PW is a mixture of dispersed oil and other hydrocarbons, dissolved inorganic and organic compounds, corrosion and treatment chemicals, heavy metals, solids, and residual production-line additives [1,2]. Produced water typically encompasses toxic pollutants, posing a tremendous environmental threat and precluding additional production activities. Taking tangible steps to treat or reuse produced water for sustainability and environmental preservation is becoming indispensable. However, PW composition and properties vary greatly based on various parameters, including the geographic field location, reservoir extraction age, reservoir geological formation, exploitation method, type of production facilities, and reservoir composition [3]. Membrane filtration has become a promising technique to meet the new environmental regulations and laws compared to earlier conventional treatment methods, such as adsorption, chemical oxidation, chemical precipitation, electrodialysis, evaporation, air flotation, gravity separation, and hydrocyclones [4]. However, membrane fouling is a significant drawback to membrane filtration due to the accumulation of oil droplets on the membrane surface and within the membrane pores. Fouling development causes concentration polarization and enhances the formation of the cake layer, which reduces the membrane performance and accelerates its deterioration [5]. Operating the membrane filtration system below the critical flux will significantly mitigate membrane fouling. Critical permeation flux is the flux at which the oil droplets' convection towards the membrane surface and their back-diffusion forces to the bulk are equal. Nevertheless, operating the membrane under the subcritical flux generates a low permeate flow [6,7].
Numerous methods have been used to control membrane fouling in crossflow filtration systems. The optimization of the operating conditions (i.e., crossflow velocity and transmembrane pressure) [8], feed turbulence promoters (i.e., reducing the concentration polarization) [9], high shear crossflow, and membrane module-surface hydrodynamic modification are examples of the conventional methods used to combat fouling development. Membrane cleaning methods have also been extensively used for fouling control; cleaning techniques can be categorized into chemical and physical processes. Chemical cleaning uses agents (e.g., alkalines and acids) to regenerate the membrane surface and clean the irreversible membrane fouling [10]. Physical cleaning uses mechanical and hydrodynamic forces to help displace and transport the oil droplets by the crossflow mainstream [11]. Mechanical cleaning techniques include ultrasonic and vibration [12]. Hydrodynamic cleaning consists of backflush pulsing using forward washing, pneumatic (air flushing, airlifting, and air bubbling), and feed pressurizing/depressurizing techniques [13].
The proposed Back-Pulsatile (BP) process (or Back-Pulsing) is a promising in situ physical antifouling technique for ceramic membrane fouling mitigation during crossflow filtration. Pulses are cyclically generated when the transmembrane pressure reaches a critical level under a constant permeate-flux filtration or at a specific time interval to reverse the filtration, loosen, and displace the oil droplets present at the membrane surface and in its pores [14]. The use of pulses in crossflow involves switching the TMP typically for one second or less at a high frequency of approximately 0.1 to 2 Hz [15]. Many research groups have extensively studied back-pulsing efficiency in crossflow filtration to control and mitigate membrane fouling while maintaining high permeation flux over the filtration time. Salalahi et al. [16] evaluated the effect of feed parameters (e.g., oil content and oil composition) and membrane flux on the efficiency of back-pulsing on ceramic-membrane filtration performance. The results showed that applying the back-pulsing significantly minimizes the fouling rate and improves the membrane performance. Gao et al. [17] examined the back-pulsing parameters (amplitude, duration, and frequency) and their interactions with membrane performance. The findings demonstrated that amplitude was the most critical variable for fouling clearance and the overall net-specific flux. However, the frequency was found to be the most important variable for the membrane's net production. Borujeni et al. [18] evaluated the performance of hollow-fiber ultrafiltration membranes at various back-pulse frequency, amplitude, and duration settings. Their results reveal that the back-pulsing operation throughout a filtration process can yield a constant flux with a substantially higher permeate-solution recovery. Ramirez et al. [19] tested the rapid back-pulsing effectiveness in improving the dilute oil-in-water clay-membrane crossflow microfiltration. The experimental results showed that the rapid back-pulsing presented high efficacy in enhancing the permeate flux approximately 25 times under an optimum duration and frequency. However, they recommended coupling the back-pulsing with other techniques (e.g., direct firing of the ceramic membrane) to maintain a higher membrane performance. Mores et al. [20] also investigated the effect of rapid back-pulsing in crossflow microfiltration on membrane performance. Direct visual observation has been accomplished to determine the development of the membrane surface fouling. The results revealed that the membrane was cleaned efficiently under a longer back-pulsing duration and higher back-pulse pressure.
The pulsatile method for treating the membrane fouling in the crossflow mode has been studied under the influence of the pulsatile parameters (amplitude, duration, and frequency). However, there are significant scientific gaps in the design studies of combining crossflow and pulsatile parameters to mitigate membrane fouling as the pulsatile process is implemented in crossflow filtration. A solid experimental design methodology should be implemented to shed light on the impact of crossflow, transmembrane pressure, back-pulse frequency, amplitude, and duration on the membrane performance.
In this study, we have assessed the crossflow filtration with pulsatile cleaning for the treatment of synthetically produced water prepared with oil from the Bakken reservoir (Southern Saskatchewan, Canada) using a ceramic membrane with an active layer of titania and a supporting layer of zirconia. Experiments were performed to measure the impact of pulsatile process and crossflow parameters (crossflow, transmembrane pressure, interval/frequency, and duration) on the membrane's overall net flux enhancement and membrane fouling control for a feed oil content of 200 ppm, which is the concentration of the Bakken reservoir-produced water after initial treatment.
A Box-Behnken response surface approach [21][22][23] was used to identify factors affecting the membrane fouling and to optimize the overall permeate flux. Three independent factors were carefully studied: the crossflow velocity (0.5 to 1.5 m/s) [24]; transmembrane pressure (0.5 to 1.5 bar) [25]; and pulsatile time interval (60 to 120 s) at a constant pulse duration (0.8 s). The overall permeate volume was collected, and its oil content and turbidity were analyzed for membrane rejection performance. The independent variable interactions and their impact on the membrane fouling were investigated to quantify their effects on the permeate flux and overall permeate volume. Finally, the pulsatile reversal flow mechanism was discussed to improve the understanding of this physical cleaning technique. The pulsatile flow efficiency was also examined to investigate the membrane performance with versus without the pulsatile process.

Experimental Methodology
In this study, a LabBrain filtration unit manufactured by LiqTech International ( Figure S1 (Supplementary Materials)) was used to perform all the experiments in the crossflow mode. The experiments were operated in batch mode, where the feed/concentrate is continuously recirculated through the ceramic membrane for 2 h. The permeate was collected, and reverse osmosis water was added to the feed to maintain the feed concentration constant. The filtration unit setup automatically logged all operational conditions such as transmembrane pressure (TMP), crossflow velocity (CFV), retentate temperature, the valve opening percentages, permeate flow rate, retentate flow rate, and feed flow rate every 3 s. For each experiment, a synthetic feed of 24 L of produced water was prepared and immediately used. The feed and permeate were characterized by measuring their oil content. Additional feed characterization was performed to measure the feed mean droplet size, zeta potential, chemical oxygen demand (COD), turbidity (TNU), pH, density, and viscosity. The ceramic membrane morphology, pore size, molecular weight cut-off (MWCO), retention capability, permeability, ceramic membrane geometry/dimensions, and thermal/chemical resistance were all reported. After each experiment, the filtration unit was meticulously washed/drained by circulating RO water to remove any impurities/ pollutants. At the end of the experiment, the ceramic membrane was cleaned using RO water (at 50 • C), Sodium hydroxide (at 50 • C), and Phosphoric acid (at 85 • C), respectively, reaching a flux recovery of 99%. The purities of the chemicals used are reported in Table 1. by an ultraviolet (UV) water purification system (EMD Millipore, Burlington, MA, USA, 2012) to prepare ultrapure deionized water (DI, <5 ppb TOC and <0.1 colony-forming units of a microorganism/mL). Sodium dodecyl sulfate (SDS, 99 wt.% pure) was purchased from Sigma-Aldrich (St. Louis, MO, USA) and used for stability in feed preparation. The 7channel ceramic membrane was purchased from Tami industries and cut into 25 × 305 mm pieces. Additional chemicals (Table 1) were purchased and used as received for the ceramic membrane cleaning or oil/solvent extraction.

Feed Synthesis and Characterization
The synthetic-produced water feed of 200 ppm ( Figures S2 and S3) was prepared by adding 4.5 mL of Bakken oil in 2 L batch of reverse osmosis water and 0.3 mM Sodium dodecyl sulfate (SDS) for feed stability. The light oil properties were measured using an Anton Paar 5000 DSA digital densitometer and a Brookfield viscometer DV-II + Pro at 22.5 • C, giving a density of 0.8872 g/cc (±5 × 10 −5 g/cc accuracy) and a viscosity of 5.23 cp (±1.0% accuracy). The 2 L oily water emulsion was mixed for 2 min at 18,000 rpm using a commercial blender to homogenize its composition and ensure feed stability. This process was repeated to prepare a total volume of 24 L (12 batches) for each experiment. Additional feed properties (Table 2) were measured using the instruments listed in Table 3.

Ceramic Membrane Characterization
This study used an ultrafiltration ceramic membrane (zirconia/titania) with seven channels (Figure 1) for the produced water treatment using the LabBrain filtration unit. Table 4 provides the essential membrane properties and specifications. The ceramic membrane active filtration and cross-sectional area were found to be 0.04186 ± 0.006 and 0.001172 ± 0.006 m 2 , respectively. The synthetic-produced water emulsion meant that the oil droplet size was measured with a Zetasizer to be around 5.2 ± 0.1 µm. In addition, the oil and water droplet adsorptions on the ceramic membrane surface were performed to measure the surface wettability responses. Figures S4 and S5 illustrate the membrane surface oleophobicity and super-hydrophilicity characteristics with contact angles of 135 and 35 • , respectively.

Ceramic Membrane Characterization
This study used an ultrafiltration ceramic membrane (zirconia/titania) with seven channels (Figure 1) for the produced water treatment using the LabBrain filtration unit. Table 4 provides the essential membrane properties and specifications. The ceramic membrane active filtration and cross-sectional area were found to be 0.04186 ± 0.006 and 0.001172 ± 0.006 m 2 , respectively. The synthetic-produced water emulsion meant that the oil droplet size was measured with a Zetasizer to be around 5.2 ± 0.1 μm. In addition, the oil and water droplet adsorptions on the ceramic membrane surface were performed to measure the surface wettability responses. Figures S4 and S5 illustrate the membrane surface oleophobicity and super-hydrophilicity characteristics with contact angles of 135 and 35°, respectively.

Crossflow Filtration System Description
The LabBrain membrane crossflow filtration system ( Figure S6) manufactured by LiqTech International (Hobro, Denmark) was used to run all the experiments under batch mode conditions. Figure S6 shows a drawing of the LabBrain Proportional and Integral (P&I) control-loop diagram system. The LabBrain unit was outfitted with a membrane housing with a ceramic element measuring 25 × 305 mm ± 1 mm. The filtration unit valves were controlled by a PLC (Siemens 6ES7 214-1AE30-OXBO, Siemens, Munich, Germany) and fed by a MotoMaster air compressor of 2.5 US gallons capacity operating at a pressure of 6 bar. The feed circulation was performed by a pump (Grundfos CRN 3-6, Grundfos, Bjerringbro, Denmark) with a 5 m 3 /h capacity at 2.5 bar. When the Data log was enabled, all parameters relating to the valves' actual openings, pressure and flow transmitters, pump speed, transmembrane pressure, flux, and production were recorded in the unit's internal memory every 3 s.
Before the start of the experiment, the ceramic membrane was progressively immersed in deionized water and then steeped for twelve hours to remove the trapped air in the membrane and to enhance the overall permeability of water flow. The membrane element was then mounted and sealed inside the housing cell. The feed emulsion was synthesized and poured into the oily wastewater tank to run the experiment. The pulsatile permeate was then collected and quantified, and the retentate and permeate were fully returned to the feed tank. Feed concentration was stably maintained by constantly adding reverse osmosis water to the container in the same sampled permeate volume during the twohour experiment.
The LabBrain filtration unit worked at settings of 1 m/s crossflow (CFV) and 1.5 bar trans-membrane (TMP). The speed of the pump and concentrate valve were progressively adjusted. Once the pressure in the unit reached a steady level, data were recorded. At the end of the experiment, a 4 mL permeate sample was characterized by measuring its oil content and turbidity.

Box-Behnken Design
The Box-Behnken design (BBD) [26][27][28] was used to study the impact of the pulsatile process in crossflow filtration mode on the membrane fouling control and mitigation. The combination of the crossflow operating conditions and pulsatile parameters was selected to understand their influence on the membrane performance. The independent process variables were the TMP, CFV, and pulsatile parameters (BP cycle). The range of the independent factors was determined based on the configuration of the experimental setup and initial single-variable experiments [6]. The pulsatile duration was constant at 0.8 s, the BP cycle varied from 60-120 s, the TMP range of 0.5 to 1.5 bar, and CFV from 0.5 to 1.0 m/s. Table 5 represents the independent factors and their coded levels −1, 0, and 1 at high, middle, and low, respectively. Permeate flux (LMH) and permeate volume (L) were selected as responses to the B-B design ( Table 6). Seventeen experiments were randomly generated with 5 central replications points and 12 design points (Table 7). In Equation (1), N is the total number of experiments, k denotes the number of investigated factors, and C 0 is the number of central replication points:  The Design Expert program (v.12) was used to analyze the responses using a quadratic polynomial regression model (Equation (2)) to fit and optimize the experimental data: where Y k is the predicted response, namely Y 1 for membrane permeate flux and Y 2 for membrane permeate volume. x i denotes the coded variables, b 0 is a model constant, and S e is the statistical error. In the equation, b i , b ii , and b ij are the linear, quadratic, and interactional regression coefficients, respectively. The model coefficients were estimated by multiple regression analysis. The polynomial model fitness was conducted using the lack of fit and the coefficient of determination (R 2 ). The experimentally examined responses were J ni (final steady permeate flux, in LMH) based on Equation (3): where J i is the flux at any time; and J 0 is the initial flux of the membrane filtration (Table S1). Y ni is the experiment's final net permeate volume measured in liters (L) calculated using Equation (4).
The membrane rejection was calculated using Equation (5): where C p and C f are the permeate and the feed oil content, respectively, in ppm; the membrane performance and fouling control were measured in each experimental run based on the permeate flux decline, overall permeate volume (Table 7), rejection capacity, and permeate turbidity (Table 8).  Figure S7 depicts the normal filtration and pulsatile flux-decline patterns. The decline in the permeate flux can be categorized into three sections. In the first section, the flux drops from its highest value at the early filtration time before reaching a pseudo-stabilized stage. Normal filtration flux decline is more significant than pulsatile flux initial decay. Consequently, the membrane surface was rapidly covered by oil droplets if the pulsatile cycle was initiated. In normal filtration mode, the crossflow field could not evacuate and transport oil droplets driven by the permeation drag field [29]. While under the pulsatile process, the decline was not very pronounced at the early stages due to the additional periodic pulsatile-cleaning impact.

Results and Discussion
The second section showed a pseudo hydrodynamic balance between crossflow cleaning and fouling development. In normal filtration mode, at 60% of initial flux, the permeate flux started to flatten until reaching a low steady-state permeate flux of 35 LMH at the final stage. Similarly, the flux decline was less manifest in the pulsatile process than in the normal filtration. The best performance of the pulsatile process was found to be for a permeate flux of 298 LMH under a periodic pulsing of a 60-s cycle, a TMP of 1.5 bar, and a CFV of 0.75 m/s. A high pulsatile permeate flux decline of 85 LMH was observed in a pulsatile cycle of 90 s, TMP of 0.5 bar, and CFV of 0.5 m/s. Table 7 depicts that an increase in the TMP to a level of 1.5 bar enhances the permeate flux and net permeate volume (Runs 1, 4, and Runs 2, 3). In addition, a high transmembrane pressure of 1.5 bar and a pulsatile reversal-flow level of 90 or 60-s lead to a high membrane performance (Runs 1, 3, 8, and 11). The central point runs (Runs 9, 12, 14, 15, and 17) at a pulsatile cycle of 90-s show a better membrane performance than all pulsatile cycles of 120 s configurations.
At the high pulsatile flow sequence of 60-s, the increase in CFV enhanced the membrane performance (Runs 5 and 6). In contrast, the increased CFV did not improve the membrane Membranes 2022, 12, 1198 9 of 21 performance at the low pulsatile frequency of 120 s. As a result, we conclude that a combination of the CFV and pulsatile flow plays a crucial role in mitigating membrane fouling. This cleaning process prevents the oil droplets, dragged by the permeation flux, from aggregating at the membrane surface and displacing them by the crossflow field, which maintains the cleanliness of the membrane surface and combats the membrane fouling.
In conclusion, the membrane steady-state permeate flux was enhanced, reaching a value of 298 LMH using the pulsatile cycle compared to a normal filtration value of 35 LMH. This can be explained by the contribution of the pulsatile process as an additional implemented physical cleaning process to combat membrane fouling and enhance membrane performance (Table 7).
Rejection improved from 92% to 99% by applying the pulsatile process. The impact of the pulsatile cycle for each run demonstrated that oil rejection for all configurations was better than in the normal filtration mode. The reduction in permeate turbidity was also enhanced from 4.58 NTU to 0.32 NTU when using this physical cleaning process (Table 8).

Box-Behnken Model
The Box-Behnken design method was used to model the 17 runs for both the experimental and predicted results, as reported in Table 3. Second-order regression polynomial equations were generated to elucidate the liaison between the membrane permeate flux/permeate volume and design-independent factors. The predictive response models were defined in coded and actual factors (Equations (6)-(9)) [30]. A list of coefficients are presented in Tables 9 and 10. The coefficient correlation value represented the corresponding factor unit change when all other factors were maintained as constant. In contrast, the intercept estimated the overall mean response of all 17 runs [31].
The obtained coefficients in Equations (6) and (7) were used for response prediction regarding flux and permeate volume, respectively. A response value was estimated for each combination of coded factor levels (1, 0, −1). In contrast, the polynomial factor coefficients identified the amplitude factor impact on the corresponding response. For the actual response quadratic regression Equations (8) and (9), the response value was predicted using the actual level unit values (TMP: 1.5,1, 0.5; CFV: 1, 0.75, 0.5; and Pulsatile cycle: 120, 90, 60), for each factor.
The full response quadratic regression equations in terms of actual factors before reduction are: The predictive model graphs for the membrane flux and permeate volume were illustrated in Figure 2a,b, respectively. The positive linear regression between all sets of experimental data for flux and permeate volume and the responses of the second-order regression polynomial equations showed that the prediction models approximated the actual data very well. The fitting statistical results between the actual and predicted values are shown in Table 11. The high values of the coefficients of determination (R 2 ) of 0.995 and 0.996 for the flux and permeate, respectively, prove a strong correlation between the experimental and predicted data. The predictive R 2 of 0.952 and 0.946 show excellent agreement with the adjusted R 2 of 0.989 and 0.990 for the flux and permeate models, respectively. The adjusted R 2 of 0.989 and 0.990 shows a high predictive capability of the two regression models for the flux and permeate volume, respectively [32]. * Adeq precision is used to measure the signal-to-noise ratio. A value higher than four means adequate model discrimination.

Analysis of Variance (ANOVA) for the Quadratic Model
The regression analysis results and the interaction between the actual factors were performed using ANOVA for membrane flux and permeate volume (Tables 12 and 13). The ANOVA of the two quadratic polynomial models was highly significant in fitting the regression models. The large model F-values of 159.63 and 179.28 indicated that there was only a 0.01% chance (p-value < 0.0001) that a higher F-value might occur due to noise for both models dealing with the membrane flux and permeate volume, respectively. The P-values of the two models were 0.0001, a very small value compared to the significance level (α = 0.05); this implies that the quadratic regression model equations showed a high goodness level of fit of the experimental response data.  For the membrane permeate-flux goodness of fit model, at a 95% confidence level, the terms with p-values less than 0.05 were significant model terms. In this scenario, A, B, AC, BC, A 2 , B 2 , and C 2 were considered. Similarly, the permeate volume model's significant terms were A, B, C, BC, A 2 , B 2 , and C 2 [33].
The F-value was used to rank the significant regression factors. The significant terms are associated with higher F-values. The following are the rankings of significant regression terms for the two responses models based on F-values: Membrane permeate-volume model: After the elimination of the non-significant interaction terms (p-value > 0.05), the final simplified-response quadratic regression models at a 95% confidence level were as follows: The full simplified response quadratic regression equations in terms of coded factors after reduction are: For all the experimental runs, the simplified quadratic regressions' Equations (10) and (11) in terms of coded factors, and Equations (12) and (13) in terms of actual factors were used to predict the permeate flux and net permeate volume, respectively.

Lack of Fit Value
As reported at the bottom of Tables 12 and 13, the results of lack of fit were not significant for both regression models relative to the pure error. The F-value for the flux model was 1.88; this implied the model's high capability to predict the experimental data variation. Similarly, the lack of fit F-value for the permeate volume model was 4.15, indicating the model's ability to estimate the considered empirical data variance [34]. In addition, Tables 12 and 13 illustrated that the relative p-values for both responses, flux and permeate volume, were 0.2747 and 0.1015, respectively, which were higher than 0.05. This meant that the models correlated with the experimental data very well, and there was a non-significant lack of fit [35].

Residual Analysis
The average probability plot residuals can also assess the model quality of fit. The plots in Figure 3 provide valuable statistical evidence. The deviation between the observed and predicted values is referred to as the residual. The normal plot of residuals is an excellent indicator of the correctness of the overall model represented by Equations (6)- (9), and defined by the Box-Behnken experimental design. The residuals in the normal probability plot were linear, indicating that the residuals satisfied the normal distribution assumption. As a result, all predicted and experimental results were in excellent agreement, validating the significance of the quadratic models presented in Equations (10)-(13) [32,33,[35][36][37].

Factors Interaction Effect and Membrane Fouling
The impact of the operating parameters (TMP, CFV, and pulsatile cycle) on the membrane permeate flux and permeate volume are presented in Figures 4 and 5, respectively. The factors' perturbation effects and their deviation from the reference configuration point (TMP: 1 bar, CFV: 0.75 m/s, and Pulsatile cycle: 90 s) are illustrated. Figure 4 demonstrates the responses' variation as a one-factor level with changes from low to high: −1, 0, +1, while the other factors were maintained constant at their specific levels. Each factor used in this experiment had a distinct impact on the membrane permeate flux and the overall permeate volume. When the CFV and Pulsatile cycle were held at their optimum level (coded value = 0), the TMP dominated the response values. The increase in the membrane permeate flux is proportionally associated with the increase of the TMP factor from 0.5 bar (coded value = −1) to 1.5 bar (coded value = +1), as shown in Figure 4a; in comparison, the membrane permeate volume stabilized at its optimum value with the increase of the TMP from 0.5 bar (coded value = −1) to 1.5 bar (coded value = +1) (Figure 5a). One should consider that a part of the permeate volume was used continuously for cleaning to maintain the membrane surface, which ensures a higher permeate flux.

Factors Interaction Effect and Membrane Fouling
The impact of the operating parameters (TMP, CFV, and pulsatile cycle) on the membrane permeate flux and permeate volume are presented in Figures 4 and 5, respectively. The factors' perturbation effects and their deviation from the reference configuration    The TMP and Pulsatile cycle were maintained at their middle levels, of 1 bar, and 90 s, respectively. The membrane permeate flux and permeate volume values improved with the increase of CFV from 0.50 cm/s (coded value = −1) to 0.75 cm/s (coded value = 0) (Figures 4b and 5b). Similarly, the pulsatile level change from 60 s (coded value = −1) to the middle level of 90 s (coded value = 0) enhanced the membrane permeate flux and permeate volume (Figures 4c and 5c), simultaneously [38].
The elliptical contour plots represent the interaction of two factors and their impact on the response value (Figures 6 and 7) while maintaining the third variable constant at its central point value. Figure 6a shows that at the central point of the pulsatile cycle at 90 s, the interaction of TMP and CFV is not significant at 95% CI. However, the increase in the range variation of the CFV and TMP or fluctuation of the design point from the central point value of the pulsatile cycle may improve the membrane permeate flux. Figure 6b illustrates the significant interaction of the TMP and Pulsatile cycle at the central point value of the CFV at 0.75 m/s. The plot shows that when TMP is at its upper range, between 1 and 1.5 bar, and the pulsatile cycle is below its central point value, the membrane permeate-flux response value was at its optimum. In comparison, the response value is low when the TMP design configuration is between 0.5 and 1 bar. In Figure 6c, the results indicate that when TMP is maintained at its central point value of 1 bar, the impact of the CFV and the pulsatile cycle is optimum at their central point values of 0.75 m/s and 90 s, respectively. A standard error calculation of the regression was performed to quantify the accuracy of the estimates. Figures S8 and S9 represent the standard error for Figures 6 and 7, respectively. The plot of the standard error of membrane permeate flux ( Figure S8) and permeate volume ( Figure S9) as a function of the independent factors (TMP, CFV, and pulsatile cycle) is presented. The smaller the standard error, the more accurate is the prediction of the regression models. In Figures S8 and S9, the darker the shading, the greater is the standard error for a given coordinate of an independent-factor response surface.   Figure 7a,b demonstrated a lack of positive interaction, at a 95% confidence interval, between the TMP and CFV at a Pulsatile cycle of 90 s and between the TMP and Pulsatile cycle at a CFV of 0.75 m/s, respectively. This was due to the dominant impact of TMP on the filtration process. At low pressure, the permeate volume reached a higher value of 8.5 L at the pulsatile cycle of 90 s and 8.8 L at a CFV of 0.75 m/s for all values of the CFV and Pulsatile cycle, respectively. However, when the TMP was set at a higher range beyond the central point value of 1 bar, the permeate volume declined for the CFV and Pulsatile variation, respectively. This could be explained by fouling development at the membrane surface due to the dominant impact of the permeation drag force resulting from the applied pressure-driven difference [39,40].
Membrane permeate volume was studied within a range of CFV from 0.5 to 1 m/s, a Pulsatile cycle from 60 to 120 s, and at 1 bar TMP. The results show that the maximum response was achieved at the middle point values of 0.75 m/s and 90 s for the CFV and Pulsatile cycle, respectively (Figure 7c).
As a result, at 95% CI, the optimum membrane-permeate volume is attained when the cleaning mechanism (CFV, Pulsatile cycle) was in the middle operating value and TMP is below its central point value. Figure 8a-c illustrates the case study for TMP at 0.7, 0.6, and 0.5 bar, respectively.

Optimization and Model Post Analysis
As reported earlier (Figures 4 and 5), each independent factor had a specific impact on the response parameters. In contrast, the factors' interactions may enhance one response but limit the other. The optimal response values were achieved by combining the process parameters levels (TMP, CFV, and Pulsatile cycle) that simultaneously satisfy the maximization criteria for each response. As a result, the desirability function was used as a multi-criteria methodology to optimize this process, identify the optimum parameters, and maximize the response values in Figures 9 and 10, as reported in Tables S2 and S3. All the independent factors and process responses are given a significance level of three and optimized in their range for the maximum membrane-permeate flux and permeate volume, respectively. The results of the numerical optimization are presented in Table 14.
The optimized process filtration, called the best run, for predicted permeate flux and membrane-permeate volume responses were 293 LMH and 8.13 L, respectively. The maximum overall desirability obtained was 0.75 for the optimum independent actual values of 1.5 bar, 0.71 m/s, and 85 s, respectively, for the approximate coded factors' levels of (1, 0, 0). An experimental validation test of the optimal values was performed to examine the results of the predicted response ( Table 14). The results showed that when the values of each factor were set at the optimal levels, the predicted and measured permeate flux and permeate volumes were exactly 293 LMH, 8.13 L, and for a desirability factor of 3, the values were 302 LMH, 8.11 L, respectively. A post-analysis has been performed running a normal filtration with no pulsatile process at 1.5 bar and 0.71 m/s for a fair comparison. The obtained permeate flux and volume results were 32 LMH and 7.13 L, respectively. Compared to the same configuration of operating pressure and crossflow velocity of 1.5 bar and 0.71 m/s with a pulsatile cycle of 85 s, the pulsatile process improved the permeability and permeate volume by 89% and 12%, respectively [14,[41][42][43].

Conclusions
This work demonstrated that pulsatile physical cleaning in crossflow filtration systems enhanced membrane filtration and fouling mitigation for a UF ceramic membrane with MWCO of 150 kg/mol. The average feed-oil content of 198 ppm was treated to produce an effluent of less than 3.5 ppm. Membrane fouling was significant without the physical cleaning. However, the implementation of the pulsatile process controls the fouling development. The cleaning performance varied greatly depending on the pulsatile process settings.
The optimization of the pulsatile process in the crossflow filtration system was successfully obtained using the Box-Behnken design method. The established regression models were verified using the analysis of variance (ANOVA) to improve the understanding of the pulsatile mechanism and to investigate the significant impact of the pulsatile cleaning technique on the membrane filtration performance and fouling control. The TMP displays the most substantial effect where a significant interaction occurred between the CFV and pulsatile time for an optimum permeability and permeate volume. In contrast, the interaction effect between the TMP and the pulsatile cycle improved the permeate flux with negligible effect on the net permeate volume due to the filtrate loss caused by the process of flow reversal in the cleaning. The optimized filtration process was achieved at a TMP of 1.5 bar and around the cleaning parameters' central points, namely a CFV of 0.71 m/s and a pulsatile time of 85 s, respectively.
Moreover, the predicted values for optimal efficiency were experimentally tested and validated for the ceramic membrane, obtaining a steady-state permeability with a flux of 302 LMH and a permeate volume of 8.11 L. In conclusion, implementing the pulsatile process improved the membrane performance and fouling alleviation in a produced water crossflow-filtration system.  Figure S1. LabBrain membrane filtration system; Figure S2. From left to right, Bakken oil and produced water feed with 200 ppm oil; Figure S3. Synthetized produced water (oil, water, and SDS surfactant); Figure S4. Contact angles of a Bakken oil droplet at the ceramic membrane surface; Figure S5. Contact angles of a water droplet at the ceramic membrane surface; Figure S6. LabBrain filtration unit P&I diagram system; Figure S7. Membrane permeate flux as a function of time for normal filtration (No pulsatile) (a) and 40 Pulsatile cycle at 60 s (b-e), 90 s (f-j), and 120 s (k-n), respectively; Figure S8. Response surface plot of membrane permeate-flux standard error as a function of TMP, CFV, and Pulsatile cycle; Figure S9. Response surface plot of membrane permeate volume standard error as a function of the TMP, CFV, and Pulsatile cycle.
Author Contributions: M.E.: writing-original draft, writing-review and editing, data curation, investigation, methodology, and conceptualization; A.H.: formal analysis, review, editing, validation, supervision, laboratory equipment, and chemicals; A.S.: review, editing, validation, and supervision. All authors have read and agreed to the published version of the manuscript.