*Article* **Through Di**ff**usion Measurements of Molecules to a Numerical Model for Protein Crystallization in Viscous Polyethylene Glycol Solution**

**Hiroaki Tanaka <sup>1</sup> , Rei Utata <sup>2</sup> , Keiko Tsuganezawa 2,3, Sachiko Takahashi <sup>1</sup> and Akiko Tanaka 1,2,3,\***


**Abstract:** Protein crystallography has become a popular method for biochemists, but obtaining high-quality protein crystals for precise structural analysis and larger ones for neutron analysis requires further technical progress. Many studies have noted the importance of solvent viscosity for the probability of crystal nucleation and for mass transportation; therefore, in this paper, we have reported on experimental results and simulation studies regarding the use of viscous polyethylene glycol (PEG) solvents for protein crystals. We investigated the diffusion rates of proteins, peptides, and small molecules in viscous PEG solvents using fluorescence correlation spectroscopy. In highmolecular-weight PEG solutions (molecular weights: 10,000 and 20,000), solute diffusion showed deviations, with a faster diffusion than that estimated by the Stokes–Einstein equation. We showed that the extent of the deviation depends on the difference between the molecular sizes of the solute and PEG solvent, and succeeded in creating equations to predict diffusion coefficients in viscous PEG solutions. Using these equations, we have developed a new numerical model of 1D diffusion processes of proteins and precipitants in a counter-diffusion chamber during crystallization processes. Examples of the application of anomalous diffusion in counter-diffusion crystallization are shown by the growth of lysozyme crystals.

**Keywords:** protein crystallization; nucleation; viscosity; diffusion; PEG; FCS; counter-diffusion; self-searching; crystallization scenario

#### **1. Introduction**

In protein-crystallization studies, polyethylene glycols (PEGs) are one of the main types of precipitants added to protein solutions to reduce their solubility and produce the supersaturation conditions required for the nucleation and growth of crystals. Many protein structures have been determined by growing protein crystals in solutions containing PEGs of various molecular weights (MW) (200–20,000 g/mol) and various concentrations (up to 40 *w*/*v*%) as precipitants [1–3]. McPherson and Gavira reported that PEGs with an MW ranging from 2000 to 8000 are the most useful precipitants, and most protein crystals were grown in solutions containing from 4 to 18% PEG for crystallization [1]. These PEG solutions are viscous, and thus may affect the mass transportation, nucleation and growth processes of the protein crystallization. However, only a few papers have experimentally studied the effect of viscosity on protein crystallization [4,5]. Therefore, we precisely measured diffusion rates in viscous solvents and found the advantages of crystallization conditions with pre-mixing viscous high-MW PEGs in the protein solution and diffusing salt as a second component of the precipitant for counter-diffusion crystallization.

**Citation:** Tanaka, H.; Utata, R.; Tsuganezawa, K.; Takahashi, S.; Tanaka, A. Through Diffusion Measurements of Molecules to a Numerical Model for Protein Crystallization in Viscous Polyethylene Glycol Solution. *Crystals* **2022**, *12*, 881. https:// doi.org/10.3390/cryst12070881

Academic Editor: Abel Moreno

Received: 30 May 2022 Accepted: 20 June 2022 Published: 21 June 2022

**Publisher's Note:** MDPI stays neutral with regard to jurisdictional claims in published maps and institutional affiliations.

**Copyright:** © 2022 by the authors. Licensee MDPI, Basel, Switzerland. This article is an open access article distributed under the terms and conditions of the Creative Commons Attribution (CC BY) license (https:// creativecommons.org/licenses/by/ 4.0/).

Crystallographic technologies are rapidly changing due to the introduction of automation and high-throughput approaches [6,7]; however, it is necessary to reach a supersaturation state for protein crystallization, so that hundreds of crystallization conditions are often tested to allow for a single protein to acquire diffraction-quality crystals. However, some important proteins are difficult to purify in large quantities. Thus, automated microfluidic systems that produce crystals using a counter-diffusion technique have been developed to save the protein sample and obtain good crystals [8]. Counter-diffusion (also known as liquid–liquid diffusion) is a common crystallization method, in which the protein and precipitants are loaded on opposite sides of a tubing chamber and gradually mixed through their diffusion [3,9,10]. In most cases, the protein is loaded into the capillary, the end of which is sealed with a gel plug to keep the protein from flowing out and allow the precipitant to diffuse in from the reservoir. When the precipitant concentration is sufficiently high and the protein chamber is sufficiently long, precipitation and crystallization of the protein occur along the tubing chamber as a result of self-searching for the optimal crystallization scenario [11]. However, in actual experiments, reagent concentrations and the lengths of the tubing chamber have limitations. Therefore, creating an equation that predicts the diffusive mass transport of proteins and precipitants in the system is crucial to optimize the crystallization conditions of counter-diffusion systems [12,13].

The translational diffusion of molecules in solution is described by the well-known Stokes–Einstein (SE) equation:

$$D = \,^k\_B T \Big/ \text{6}\pi\eta r \tag{1}$$

where *k<sup>B</sup>* is the Boltzmann constant. The diffusion coefficient (*D*) depends on the hydrodynamic radius of the molecule approximated as a sphere with radius *r*, temperature *T*, and solvent viscosity η. According to the SE equation, the diffusion rates of molecules are inversely proportional to the viscosity of the solvent, although, in some macromolecular crowding solvents, such as high-MW PEGs, the measured diffusion coefficient is larger than a value calculated using the macroviscosity of the solution. The phenomena are known as anomalous diffusion [14,15]. We proposed empirical equations to estimate the *D* value of a protein in a PEG solution from the protein, PEG MWs and concentration of the PEG [16]. Later, Holyst et al. showed that, when considering the viscosity around the molecules (nanoviscosity), the SE equation is correct, even in viscous PEG solutions [17]. Thus, we planned an experimental study of the anomalous effects of PEG solvents on the diffusion rates of small molecules, peptides, and proteins to predict *D* values in viscous PEG solvents.

We used fluorescence correlation spectroscopy (FCS) to analyze the diffusion phenomenon. FCS is an effective experimental method to observe the translational diffusion of a molecule under various solvent conditions. The original concept of FCS is to detect and analyze the spontaneous fluctuations in the fluorescence emissions of several molecules in a small detection volume caused by thermodynamic fluctuation [18]. The translational diffusion time (τ*D*), in which a molecule remains in the focal volume depending on its hydrodynamic radius, and is proportional to the cubic root of its molecular mass for a spherical particle. Experimentally, τ*<sup>D</sup>* is determined by calculating the autocorrelation function from the fluorescence-intensity fluctuations caused by fluorescent molecules diffusing in and out of the detection volume.

According to Holyst et al. [17], the SE equation can be rewritten as Equation (2). If the diffusion constant of the solute in water is assumed to be *D*<sup>0</sup> and the viscosity of water η0, this can be written as:

$$D = \,^{k\_{\rm B}} \Big/ \not{\pi} \pi \eta\_{\rm nano} \tau = \,^{D\_0 \times \eta\_0} \Big/ \not{\eta}\_{\rm nano} \tag{2}$$

where *D* and η*nano* are the diffusion constant of the solute and the nanoviscosity around the solute molecule in the solvent. The τ*<sup>D</sup>* value of a solute is inversely proportional to *D*; therefore, Equation (2) shows the following relation:

$$
\tau\_{\rm D} = \frac{k\_1}{\rho} \Big/\_{\rm D} = \frac{k\_1 \eta\_{\rm namo}}{} \Big/\_{\rm D\_0 \eta\_0} = \frac{k\_1 \eta\_{\rm nono} \tau\_{\rm D0}}{} \Big/\_{\eta\_0} = k\_2 \tau\_{\rm D0} \times \eta\_{\rm nono} \tag{3}
$$

where *k*<sup>1</sup> and *k*<sup>2</sup> are specific constants.

We measured the τ*<sup>D</sup>* values of small molecules, peptides, and proteins in various viscous PEG solvents in this paper. In viscous, high-MW-PEG solutions, the diffusion rates of small molecules are selectively susceptible to anomalous diffusion and not inversely proportional to the viscosity of the solution (macroviscosity). Our previous model [16] is greatly improved in this paper with the use of the newly measured results.

The substantial difference between the diffusions of a small molecule and a protein in PEG solution indicated the following benefits to counter-diffusion crystallization in viscous high-MW-PEG solutions: the delays in the diffusion rates of the various precipitants, such as salts and organic solvents, will be smaller than previously expected, and, at the same time, the protein leakage from the crystallization chamber will be suppressed. We have developed a new numerical model of the 1D diffusion processes of the proteins and precipitants in the counter-diffusion chamber, based on the diffusion times of molecules newly measured in the PEG solutions. The new calculation allows us to point out the advantage of using high-MW PEGs in the counter-diffusion crystallization. We have also demonstrated the advantage with an example of lysozyme crystallization.

#### **2. Materials and Methods**

#### *2.1. Materials*

PEG1000, PEG3350, PEG6000, PEG10,000, and PEG20,000 were purchased from Hampton Research Chemicals (Aliso Viejo, CA, USA). The fluorescent molecules, TAMRA (MW = 528.0) and ALEXA647 (MW = 1155.1), were purchased from Olympus (Tokyo, Japan). The TAMRA-labeled synthetic peptides, SHP-1 binding peptide (10-mer TAMRA: ITpYSLLKGGK-TAMRA, MW = 1572.6) and p53 N-terminal peptide (16-mer TAMRA: SQETFSDLWKLLPEN-K-TAMRA, MW = 2346.6), were purchased from Toray Research Center, Inc. (Tokyo, Japan). Alexa Fluor 647-labelled goat anti-human immunoglobulin G (IgG, MW = 150,000) was purchased from Molecular Probes (Eugene, OR, USA). Egg-white lysozyme was purchased from FUJIFILM (Tokyo, Japan).

#### *2.2. Preparation of the Protein Samples*

Human-bromodomain-containing protein 2 (Brd2) (accession: P25440.2, 74aa-194aa, MW = 15,000), and mouse secernin-1 (accession: Q9CZC8.1, 1aa-414aa, MW = 46,300) were synthesized by the Escherichia coli. cell-free protein-synthesis method [19]. The coding cDNA fragments were amplified by polymerase chain reaction, to add a T7 promoter and histidine affinity tag-encoding sequence to the 5<sup>0</sup> -region, and a T7 terminator sequence to the 3<sup>0</sup> -region. These subclones were ligated into the pCR2.1-TOPO vector (Invitrogen, Carlsbad, CA, USA) for cell-free protein expression [19–21].

The proteins produced by the *E. coli* cell-free synthesis system with an N-terminal histidine affinity tag with a TEV cleavage site were purified by a HisTrap HP column (Cytiva, Marlborough, MA, USA). The histidine affinity tag was then removed by incubation with TEV protease at 277 K overnight, and the proteins were further purified by a HisTrap HP column and concentrated [22]. Thus, the final protein samples contained additional amino acid sequences derived from the expression vectors, that is, GSSGSSG for the N terminus and SGPSSG for the C terminus. The molecular weights of the final protein samples were 15,388.5 for Brd2 and 47,317.1 for secernin-1.

#### *2.3. Protein Labeling*

The amine side chains of the proteins were randomly labeled with Alexa Fluor 647 mono-functional succinimidyl ester (Molecular Probes, Eugene, OR, USA). The Alexa Fluor 647 fluorophore was covalently attached to the protein by the conjugation protocol of the manufacturer. In brief, 13.7 mM protein was dissolved in 0.1 M Na2CO<sup>3</sup> (pH 8.3) containing the fluorescent dye, and the mixture was stirred for 1 h at room temperature. The free dye molecules were then removed using BioGel P-30 Fine size-exclusion purification resin (Bio-Rad, Hercules, CA, USA). After the labeling step, absorption spectra indicated

that 1.26 mole Alexa Fluor was contained in 1 mole Brd2 molecules, and 1.12 mole in 1 mole secernin-1 molecules.

#### *2.4. FCS Measurements*

The FCS measurements of the solution samples were performed with an MF20 singlemolecule fluorescence detection system (Olympus, Tokyo, Japan) using the on-board 543 nm (for TAMRA) or 633 nm (for Alexa Fluor) helium–neon laser at a laser power of 100 µW for excitation [23,24]. For convenience, the experiments were performed in 384-well glass-bottom plates using a sample volume of 30 µL.

The fluorescent samples and each PEG solution (*w*/*v*) were mixed in FCS buffer containing 50 mM Tris/HCl (pH 8.0) and 0.05% Tween 20. Tween 20 was added to suppress glass–surface interactions. The final concentrations of the fluorescent molecules used for the measurements were adjusted to a concentration of 1 nM. All of the FCS measurements were performed in duplicate. The measurement data were obtained with a data acquisition time of 10 s per measurement, and the measurements were performed five times per sample at 296 K. For machine performance verification and normalization of the obtained results, the τ<sup>D</sup> values of the standard fluorescent dyes (1 nM TAMRA or 1 nM Alexa Fluor 647) were determined at each measurement. FCS data analysis was performed with the MF20 software package (Olympus). The error bars shown in the graphs represent the standard deviations of five measurements. The Brd2 and secernin-1 samples were purified after fluorescent labelling, although the samples still contained the free-labeling reagent (29% for the Brd2 sample, and 16% for the secernin-1 sample). Therefore, the τ<sup>D</sup> values were analyzed by the two-component analysis software of the instrument (MF20) and decomposed into two τ<sup>D</sup> values owing to the protein and contaminating free-labeling reagent (Alexa647). For the analysis, the τ<sup>D</sup> value of Alexa647 in each PEG solvent was measured and used.

#### *2.5. Viscosity Measurements*

Each PEG solution (*w*/*v*%) containing the FCS buffer was prepared by gentle stirring. The viscosities of the PEG solutions were measured at 296 K by a SV10 sine-wave vibroviscometer (A&D Company, Ltd., Tokyo, Japan) calibrated with JS10 and JS100 calibration solutions (Nippon Grease Co., Ltd., Yokohama, Japan).

#### *2.6. Crystallization*

The gel-tube counter-diffusion methods [13,25] were used, as previously mentioned. A capillary (0.47 mm bore) with a gel tube was filled with 7 µL of a 20% PEG solution (PEG 4000, 10,000, or 20,000) containing 20 mg/mL lysozyme and 50 mM acetate buffer (pH 4.5). The capillaries were vertically placed in a reservoir solution containing 600 mM NaCl, 20% PEG, and 50 mM acetate buffer (pH 4.5) for 30 days at 293 K. Gel tubes with a 1 mm bore were presoaked in each reservoir solution for more than 1 week before sample loading.

#### *2.7. X-ray Di*ff*raction Experiment*

The crystals grown in the capillary were picked out and immersed in cryo-protectant solution, including 600 mM NaCl, 40% PEG4000, and 50 mM acetate buffer (pH 4.5) using a cryoloop. They were then flash-frozen with liquid nitrogen.

Data collection was performed using synchrotron radiation at Diamond Light Source beamline i04 equipped with an Eiger2 XE 16M pixel detector. All of the datasets were integrated and scaled using the programs iMosflm [26] and Aimless [27], as implemented in the CCP4 program package [28].

#### **3. Results and Discussion**

#### *3.1. Macroviscosity (*η*macro) of PEG Solutions*

The PEG polymers with various MW (1000–20,000) were dissolved in an FCS buffer at 1–20% (*w*/*v*), and their viscosities were measured (Table 1). The high-MW PEGs (PEG10,000 and PEG20,000) strongly increased the viscosities of the solutions and the lower-MW PEGs (PEG1000 and PEG3350) were not effective, even at 20%. The measured densities of the solutions were 1.01 g/cm<sup>3</sup> (1, 2 and 5%), 1.02 g/cm<sup>3</sup> (10%), 1.03 g/cm<sup>3</sup> (15%), and 1.04 g/cm<sup>3</sup> (20%).


**Table 1.** Measured macroviscosities of the PEG solutions in the FCS buffer.

\* The PEG solutions contained 50 mM Tris-HCl (pH 8.0), and 0.05% Tween 20 for stable FCS measurements.

#### *3.2. Di*ff*usion Times (*τ*D) of Compounds, Peptides, and Proteins in Viscous PEG Solutions*

The effects of the various PEG solvents (PEG1000, PEG3350, PEG6000, PEG10,000, and PEG20,000 solvents) on the diffusion rates of solute molecules were measured by FCS. We measured the diffusion rates of seven solute molecules: the fluorescence dye TAMRA (Mw = 528.0), Alexa647 (1155.1), fluorescence-labeled peptide 10-mer TAMRA (1572.6), and 16-mer TAMRA (2346.6) (Figure 1), fluorescent-labeled proteins, Brd2 (16,500), secernin-1 (48,500), and IgG (150,000) (Figure 2). The Brd2 and secernin-1 samples were purified after fluorescent labelling, although the samples still contained the free-labeling reagent (29% for the Brd2 sample, and 16% for the secernin-1 sample). Therefore, the τ*<sup>D</sup>* values were analyzed by the two-component analysis software of the instrument (MF20) and decomposed into two τ*<sup>D</sup>* values owing to the protein and contaminating free-labeling reagent (Alexa647). For the analysis, the τ*<sup>D</sup>* value of Alexa647 in each PEG solvent was measured and used. The τ*<sup>D</sup>* values of each solute in the various PEG solutions were normalized by their standard τ*D*value (τ*D*0), determined in FCS buffer solution without PEG. The normalized τ*<sup>D</sup>* values are plotted as a function of the macroviscosity of the PEG solution in Figures 1 and 2.

Figures 1 and 2 show that, in the PEG1000 solution, the normalized τ*<sup>D</sup>* values of the seven molecules almost linearly increased with increasing solvent viscosity, which mostly followed the SE equation (Equation (3)). However, the slopes of the regression curves decreased with the increasing MW of the PEGs: PEG1000 (purple) > PEG3350 (blue) > PEG6000 (green) > PEG10,000 (orange) > PEG20,000 (red). In high-MW-PEG solutions, the normalized diffusion mobility was strongly nonlinear and showed considerable deviation to faster diffusion than the SE behavior. The extent of the deviation of normalized τ*<sup>D</sup>* values was more significant when the solute MW was smaller: 16-mer TAMRA (2346.6) < 10-mer TAMRA (1572.6) < Alexa647 (1155.1) < TAMRA (528.0). No clear difference was observed with Brd2 (16.5 k), secernin-1 (48.5 k), and IgG (150 k). These results showed that the extent of diffusion mobility deviations depends on both the MWs of the solute (*Msol*) and the solvent PEG (*Mpeg*). When the *Msol* is smaller and the *Mpeg* is larger, the measured τ*<sup>D</sup>* becomes smaller than expected by the SE equation. Note that the normalized τ*<sup>D</sup>* of the small molecular compounds did not exceed 6-fold, even in 20% PEG20,000 solution with a viscosity of 52.3 mPa·s.

Figure 3 shows the correlations between the normalized solute τ*<sup>D</sup>* and the solvent viscosity, which are shown in Figures 1 and 2, in a double logarithmic plot. The correlations were linearly regressed by the least-squares method and fitted well as the following:

$$
\mathcal{L}\log(\tau\_D/\tau\_{D0}) = \mathcal{L}\log A + \alpha \mathcal{L}\log \eta\_{macro}
$$

**Figure 1.** Measured and normalized τ*<sup>D</sup>* values of TAMRA, Alexa647, 10-mer TAMRA, and 16-mer TAMRA in different PEG solvents. The τ*<sup>D</sup>* values were measured and normalized by that measured in FCS buffer solution without PEG. The measurements (*n* = 5) were performed in duplicate. The purple, blue, green, orange, and red points are the results in the PEG1000, PEG3350, PEG6000, PEG10,000, and PEG20,000 solutions, respectively. The data were analyzed and plotted using GraphPad Prism version 8.4.2 for Windows (GraphPad Software, San Diego, CA, USA).

**Figure 2.** Normalized τ*<sup>D</sup>* values of Brd2, secernin-1, and IgG in different PEG solvents. The measurements (*n* = 5) were performed in duplicate. The measured τ*<sup>D</sup>* values were analyzed by the two-component analysis software of the instrument (MF20) and decomposed into two τ*<sup>D</sup>* values owing to the protein and contaminating free-labeling reagent, then calculated τ*<sup>D</sup>* values of proteins were normalized as in Figure 1.

**Figure 3.** The measured diffusion times (τ*D*) depend on macroviscosity (η*macro*). The relationships between the normalized solute τ*D*s and solvent η*macro*s are shown with both logarithmical axes for TAMRA, 16-mer TAMRA and IgG.

By substituting τ*<sup>D</sup>* by referring to Equation (3), this equation is written as Equation (4) with an anomalous index α (α 5 1).

$$\text{Log}(\tau\_D | \tau\_{D0}) = \text{Log}(\mathbf{k}\_2 \eta\_{\text{mano}}) = \text{Log}\mathbf{k}\_2 + \text{Log}\eta\_{\text{mano}} = \text{Log}\mathbf{k}\_2 + \alpha \text{Log}\eta\_{\text{macro}} \tag{4}$$

When SE equation holds α is 1 (Equation (3)), and α would not exceed 1, because the diffusion becomes faster than the SE relation. The estimated α values are 0.9235, 0.6628, 0.5807, 0.5088, 0.4843 and 0.393 for PEG1000, PEG3350, PEG6000, PEG8000, PEG10,000 and PEG20,000, respectively, for TAMRA. 0.9553, 0.7992, 0.6959, 0.6407, 0.6062 and 0.5049 for 16-mer TAMRA and 0.8564, 0.8257, 0.7264, 0.6631, 0.6511 and 0.5699 for IgG. The values of *k*<sup>2</sup> estimated from Equation (4) are listed in Table 2, (their average value is 0.907).


**Table 2.** Obtained *k*<sup>4</sup> , *RSP* and *k*<sup>2</sup> values for the solutes.

\* The average of *k*<sup>2</sup> is 0.907.

#### *3.3. Quantitative Approximation of the Anomalous Di*ff*usion*

It seems plausible that the extent of anomalous diffusion depends on the relative size of the solute molecule and the PEG molecule. Thus, in Figure 4, the anomalous index, α, is plotted against the ratio of the *Msol* and the *Mpeg*, *Msol*/*Mpeg*, with both logarithmic axes. We found that relations of all α and the ratio can be approximated by exponential relations as:

$$\log(\alpha) = \text{Log}k\_3 + k\_4 \text{Log}\left(\frac{M\_{\text{sol}}}{M\_{\text{pg}}}\right) = \text{Log}\left(k\_3 \times \left(\frac{M\_{\text{sol}}}{M\_{\text{pg}}}\right)^{k\_4}\right)$$

$$\alpha = k\_3 \times \left(\frac{M\_{\text{sol}}}{M\_{\text{pg}}}\right)^{k\_4} \tag{5}$$

where *k*<sup>3</sup> and *k*<sup>4</sup> are specific constants.

**Figure 4.** The anomalous index α depends on *Msol*/*Mpeg*. Blue, orange, grey, yellow, light blue, green and dark blue dots correspond to the α value of TAMRA, ALEXA, 10-mer TAMRA, 16-mer TAMRA, Brd2, sesernin-1 and ALEXA-IgG obtained in various molecular weight PEGs, that is PEG1000, PEG3350, PEG6000, PEG8000, PEG10,000 and PEG20,000. Colored dotted lines show exponential approximations for the corresponding dots.

We defined *RSP* as the ratio of *Msol* and *Mpeg*, where α reaches 1:

$$1 = k\_{\mathbb{3}} \times RSP^{k\_4} \tag{6}$$

and Equation (5) can be written as:

$$\alpha = \left(\frac{M\_{\rm sol}}{RSP \times M\_{\rm peg}}\right)^{k\_4} \tag{7}$$

$$\log \alpha = k\_4 \log \left(\frac{M\_{\text{sol}}}{RSP \times M\_{\text{pg}\text{g}}}\right) = k\_4 \log \left(\frac{1}{RSP}\right) + k\_4 \log \left(\frac{M\_{\text{sol}}}{M\_{\text{pg}\text{g}}}\right)$$

Then, *k*<sup>4</sup> and *RSP* were determined by calculating the data in Figure 4 using the method of least-squares. They are shown in Table 2.

The obtained *k*<sup>4</sup> and *RSP* values are plotted against *Msol* in Figure 5 with both logarithmic axes. Both parameters can be approximated by exponentiations as follows:

$$\text{Logk}\_4 = \text{Logk}\_5 - k\_6 \text{LogM}\_{\text{sol}} \tag{8}$$

$$\text{LogRSP} = \text{Logk}\gamma + k\_8 \log \mathcal{M}\_{\text{sol}} \tag{9}$$

where *k*5, *k*6, *k*<sup>7</sup> and *k*<sup>8</sup> are specific constants. The dotted lines in Figure 5 are the linear regression of the dots, showing that *<sup>k</sup>*5, *<sup>k</sup>*6, *<sup>k</sup>*<sup>7</sup> and *<sup>k</sup>*<sup>8</sup> are 0.6804, 0.1450, 3.819 <sup>×</sup> <sup>10</sup>−<sup>4</sup> and 1.1456, respectively. Then, the anomalous index, α in the Equation (7), is as follows, using the estimated values:

$$\alpha = \left(\frac{M\_{\rm sol}}{k\_7 \times M\_{\rm sol}{}^{k\_8} \times M\_{\rm peg}}\right)^{k\_5 \times M\_{\rm sol}{}^{-k\_6}} = \left(\frac{M\_{\rm sol}}{3.819 \times 10^{-4} \times M\_{\rm sol}{}^{1.1456} \times M\_{\rm peg}}\right)^{0.6804 \times M\_{\rm sol}{}^{-0.145}}\tag{10}$$

**Figure 5.** *k*<sup>4</sup> and *RSP* values plotted against *Msol*. Both parameters can be approximated by exponentiation with *Msol*.

Finally, using Equation (10), we are able to estimate the anomalous index α, so that τ*<sup>D</sup>* of a solute molecule in an arbitrary PEG solvent can be estimated by the following equation.

$$
\pi^{\mathsf{T}\_{\mathsf{D}}} \Big/\_{\mathsf{T}\_{\mathsf{D}0}} = k\_{\mathsf{2}} \times \eta\_{\text{macro}}{}^{\alpha} \tag{11}
$$

#### *3.4. Simulation of the Di*ff*usion Processes in the Counter Di*ff*usion Chamber*

It is important for counter-diffusion crystallization to estimate the changes over time regarding the concentrations of the protein and crystallization reagents in the counterdiffusion chamber. For the simulation program, we first estimated the macroviscosities of arbitrary MW PEGs with various concentrations. We have reported empirical equations to estimate the values as Equations (1) and (2) in the previous paper [16], so that the results of Table 1 were assigned to the equations, and the parameters were refined by the least-squares methods using the Microsoft Excel solver. The equations with the refined parameters are as follows:

$$
\eta\_{macro} = 1.002 \times e^{\circ \times \text{C}\_{\text{pre}}} \tag{12}
$$

$$\gamma = 0.045761 \times \ln(Mpc) - 0.2554 \tag{13}$$

where *Cpeg* is the *w*/*v*% of the PEG solvent. From Equations (12) and (13), η*macro* of PEG solvents can be calculated from their MWs and concentrations.

The changes in protein and crystallization reagents' concentrations in each area (0.5 mm length) of the counter-diffusion chamber since the start of diffusion can be estimated from each initial concentration by solving the one-dimensional diffusion partial difference equations repetitively by Microsoft Excel macro, as follows:

$$D\_{\rm ppv}(\mathbf{x},t) = D\_{\rm ppv0} \Big/ \left(0.907 \times \eta\_{\rm macro} \Big(\mathcal{M}\_{\rm peg}, \mathcal{C}\_{\rm peg}(\mathbf{x},t)\Big)^{a(\mathcal{M}\_{\rm ppv}, \mathcal{M}\_{\rm peg})}\right) \tag{14}$$

$$D\_{\rm sol}(\mathbf{x},t) = \,^{D\_{\rm sol0}} \Big/ \Big( 0.907 \times \eta\_{\rm macro} \Big( \mathcal{M}\_{\rm peg}, \mathcal{C}\_{\rm peg}(\mathbf{x},t) \Big)^{a(M\_{\rm sol}, M\_{\rm peg})} \right) \tag{15}$$

$$\text{if } Ar(\mathbf{x} - \Delta \mathbf{x}) \le Ar(\mathbf{x}) \text{ then } \beta\_1 = \frac{Ar(\mathbf{x} - \Delta \mathbf{x})}{Ar(\mathbf{x})} \text{ else } \beta\_1 = 1$$

$$\text{if } Ar(\mathbf{x} + \Delta \mathbf{x}) \le Ar(\mathbf{x}) \text{ then } \beta\_2 = \frac{Ar(\mathbf{x} + \Delta \mathbf{x})}{Ar(\mathbf{x})} \text{ else } \beta\_2 = 1$$

$$\left(-\Delta \mathbf{x}, t\right) - C\_{\text{err}}(\mathbf{x}, t)\Big|\ge \beta\_1 + \left(C\_{\text{err}}(\mathbf{x} - \Delta \mathbf{x}, t) - C\_{\text{err}}(\mathbf{x}, t)\right) \times \beta\_2\right) \times D\_{\text{err}}(\mathbf{x}, t) \times \frac{\Delta t}{\Delta t} \quad (16)$$

$$\mathbf{C}\_{\rm pwo}(\mathbf{x},t+\Delta t) = \mathbf{C}\_{\rm pwo}(\mathbf{x},t) + \left( \left( \mathbf{C}\_{\rm pwo}(\mathbf{x}-\Delta t,t) - \mathbf{C}\_{\rm pwo}(\mathbf{x},t) \right) \times \boldsymbol{\beta}\_1 + \left( \mathbf{C}\_{\rm pwo}(\mathbf{x}-\Delta t,t) - \mathbf{C}\_{\rm pwo}(\mathbf{x},t) \right) \times \boldsymbol{\beta}\_2 \right) \times \mathbf{D}\_{\rm pwo}(\mathbf{x},t) \times \frac{\Delta t}{\Delta x^2} \tag{16}$$

$$\mathbb{C}\_{\text{ssl}}(\mathbf{x},t+\Delta t) = \mathbb{C}\_{\text{ssl}}(\mathbf{x},t) + \left( (\mathbb{C}\_{\text{ssl}}(\mathbf{x}-\Delta \mathbf{x},t) - \mathbb{C}\_{\text{ssl}}(\mathbf{x},t)) \times \mathbb{} \mathbb{I}\_{1} + (\mathbb{C}\_{\text{ssl}}(\mathbf{x}-\Delta \mathbf{x},t) - \mathbb{C}\_{\text{ssl}}(\mathbf{x},t)) \times \mathbb{I}\_{2}) \times \mathbb{D}\_{\text{ssl}}(\mathbf{x},t) \times \frac{\Delta t}{\Delta x^{2}} \tag{17}$$

$$\mathbb{C}\_{\text{pg\%}}(\mathbf{x}, t + \Delta t) = \mathbb{C}\_{\text{pg\%}}(\mathbf{x}, t) + \left( \left( \mathbb{C}\_{\text{pg\%}}(\mathbf{x} - \Delta \mathbf{x}, t) - \mathbb{C}\_{\text{sol}}(\mathbf{x}, t) \right) \times \boldsymbol{\beta}\_1 + \left( \mathbb{C}\_{\text{pg\%}}(\mathbf{x} - \Delta \mathbf{x}, t) - \mathbb{C}\_{\text{sol}}(\mathbf{x}, t) \right) \times \boldsymbol{\beta}\_2 \right) \times D\_{\text{pg\%}} \times \frac{\Delta t}{\Delta \mathbf{x}^2} \tag{18}$$

where *Dpro*0, *Dsol*<sup>0</sup> and *Dpeg* are the diffusion coefficients of the protein, other solute, and PEG in water. *Dpro*(*x*, *t*) and *Dsol*(*x*, *t*) are the diffusion coefficients of the protein and other solute at position *x* in the chamber on time *t*. They are derived from Equations (10)–(13), using *Cpro*(*x*, *t*), *Csol*(*x*, *t*) and *Cpeg*(*x*, *t*) which are the concentrations of the protein, other solute, and PEG at *x* and on *t*. *Ar*(*x*) is the cross section of the chamber at *x*. β<sup>1</sup> and β<sup>2</sup> are the cross-section factors. They are determined by the ratio of the two cross sections of adjacent regions having length of ∆*x* (in this case, 0.5 mm). The calculation step time, ∆*t* is 20 s. According to Equations (14) and (15), the diffusion coefficients (*D*) of molecules (protein, other solute) in the small area at time *t* are calculated dividing the values in the aqueous solution (*D0*) by each calculated ratio, τ*D*/τ*D0*. The τ*D*/τ*D0* of molecules was calculated according to the macroviscosity of the small areas (Equations (10) and (11)), and the value of macroviscosity is determined by concentration and MW of PEG using Equations (12) and (13). The diffusion coefficients of PEG4000, PEG20,000, NaCl, and lysozyme in water are 1.24 [29], 0.476 [29], 15.0 [30], 1.06 [31], ×10−10/m<sup>2</sup> s −1 , respectively. The counter-diffusion chamber was assumed to be a gel tube with a length of 6 mm and an inner diameter of 1 mm, and a capillary with an inner diameter of 0.5 mm connected to the gel tube, and a protein sample solution was filled with 40 mm in the capillary.

Figure 6 show the results of repetitive calculation of Equations (14)–(18). Figure 6a shows the process of a case where a capillary is filled with 24 mg/mL lysozyme and a crystallization reagent is 20% PEG4000. Figure 6a1 shows the protein concentrations in various areas of the capillary after 4, 32, 64, 128, and 192 days. During this period, the lysozyme diffused out of the capillaries, and its concentration significantly decreased (Figure 6a1). Figure 6a2 shows those of PEG4000. It took time for PEG4000 to diffuse into the capillary. Figure 6a3 plots the PEG4000 concentrations of each part of the capillary on the horizontal axis and the corresponding protein concentrations on the vertical axis. The plotted curves show that the capillary could scan the protein-PEG4000 plane over time, but there is a wide unscanned area in the upper right corner of the figure. Figure 6b shows the case of 20% PEG20,000 with the same lysozyme concentration. Since the diffusion coefficient of PEG20,000 is about 1/3 of that of PEG4000, it took more time to diffuse into the capillary (Figure 6b2). During this period, the lysozyme diffused out of the capillaries, and its concentration significantly decreased (Figure 6b1). As a result, the scanned area of the protein-PEG plane becomes narrower (Figure 6b3). These results show that the self-searching mechanism for the optimal crystallization scenario would not work well when high-MW PEG is used as a precipitant in the counter-diffusion method.

These results indicated that viscous PEG solvents are not good precipitants for counterdiffusion crystallization. On the other hand, high-MW PEG is expected to improve the quality of crystals due to its dehydration effect [32]. Thus, PEG is a reagent that we would like to try.

Regarding the nucleation, the probability is explained by the following equation [33]:

$$\frac{\partial \mathcal{N}}{\partial t} = V \times \frac{const}{\eta} \times \exp\left(-\frac{16\pi\nu^2\gamma^3}{3(kT)^3[\ln \mathcal{S}]^2}\right) \tag{19}$$

where ∂*N*/∂*t, S,* γ*,* and ν*,* are the nucleation probability, supersaturation, surface energy, and volume of the crystal. *V,* η and *const* are the volume of the solution, viscosity and the constant, which are related to the attachment kinetics of growth units. This depends on the molecular charge, the molecular volume, and the density of the solution. Therefore, it is expected that the nucleation probability will decrease in highly viscous solution.

**Figure 6.** Calculated results of the diffusion processes in the crystallization chamber with 1D simulation. The results of repetitive calculation of Equations (14)–(18) are shown. The crystallization chamber is assumed to be a 6 mm length and 1 mm bore gel tube and 0.5 mm bore glass capillary with a 40 mm length protein solution; (**a**) 24 mg/mL of lysozyme is loaded into the capillary and 20% of PEG4000 is applied as the precipitant for the reservoir solution; (**b**) 24 mg/mL of lysozyme and 20% of PEG20,000; (**c**) 24 mg/mL of lysozyme with 20% of PEG20,000 is loaded into the capillary, 20% of PEG20,000 is presoaked in the gel tube and 20% of PEG20,000 and 600 mM of NaCl is applied as the precipitant for the reservoir solution. The protein concentrations along the chamber are shown in (**a1**,**b1**,**c1**). The abscissa shows the position along the chamber from the open end of the gel tube. The concentrations of the precipitant along the chamber are shown in (**a2**,**b2**,**c2**). The concentration relations of the precipitant and the protein along the chamber are plotted in (**a3**,**b3**,**c3**). For (**a**,**b**), \_, , N, ×, Ж show the values of 4, 32, 64, 128 and 192 days after the starting of the diffusion, respectively. For (**c**), they show the values of 0.2, 1, 2, 8 and 16 days, respectively. The values are plotted every 2 mm along the chamber.

We reported the nucleation probability in Table 1 of our previous report [34]. The lysozyme was crystallized by batch method with a 50 mM acetate buffer of pH 4.5 and various amounts of PEG4000 and NaCl. The concentrations of NaCl tested were 300, 500, 800 and 1000 mM, and of PEG4000, they were 0, 5, 10, 15, 20 and 25%. The nucleation probability decreased in the higher concentration PEG solvents, almost inversely proportional to the τ*D*/τ*D*<sup>0</sup> values. The result also suggested the possibility of viscous PEG solvents to obtain fewer but larger crystals.

PEG and small molecules, such as salts, synergistically act on protein crystallization [35]. Therefore, in Figure 6c, we have simulated the following case: a protein sample is pre-mixed with 20% PEG20,000 in advance and crystallization begins upon diffusing 600 mM NaCl contained in 20% PEG20,000 solvent into the capillary. The PEG20,000 concentration was pre-uniform in the counter-diffusion chamber. The concentrations were plotted 0.2, 1, 2, 8 and 16 days after the start of diffusion. The results (Figure 6c3) showed that, after 16 days, most of the protein-NaCl plane was scanned. From these results, it was strongly suggested that crystallization using the counter-diffusion method with high-MW PEG is both possible and promising.

#### *3.5. Quality of Crystals Grown in Viscous PEG Solvents Using Counter-Di*ff*usion Systems*

To confirm the crystallization performance obtained by diffusing a small-molecule precipitant into the capillary, where protein samples were pre-mixed with viscous PEG solution in advance, lysozyme crystals were grown in PEG4000, PEG10,000, and PEG20,000 solutions using the counter-diffusion method with gel-tube parts. A lysosome solution with 20 mg protein/mL in each 20% PEG solution containing 50 mM acetate buffer (pH 4.5) was filled in the capillaries and placed in reservoir solutions containing the corresponding 20% PEG, 50 mM acetate buffer (pH 4.5), and 600 mM NaCl at 293 K. After incubation for 20 days at 293 K, crystals were observed in all the capillaries and, after incubation for 30 days, the crystals were harvested. In the X-ray diffraction experiment, three or four crystals grown in each PEG solution were used and the datasets were collected. The crystal-to-detector distance was fixed to the maximum resolution of 1.12 Å. For the crystal from the PEG20,000 solution, an additional dataset with the distance fixed to 1.02 Å was also collected from the same crystal. The dataset from the crystal with the lowest *B* factor of the Wilson plot was selected, and the statistics are listed in Table 3.

**Table 3.** Data collection and scaling statistics. The values in parentheses are for the highest resolution shells.


\* The dataset was collected with the proper crystal-to-detector distance, and observed the maximum resolution of 1.02 Å.

The parameters that are usually used to evaluate the protein-crystal quality are the maximum resolution, *R*merge, <*I*>/<σ(*I*)>, *B* factor of the Wilson plot, and mosaicity [36,37]. Among these parameters, the maximum resolution, *R*merge, and <*I*>/<σ(*I*)> depend on the experimental diffraction conditions, such as the size of the crystals and intensity of the X-ray. Conversely, the overall *B* factor obtained from the relative Wilson plot is independent of the experimental diffraction conditions [36]. The Wilson *B* factor of the crystal grown in the PEG20,000 solution was the smallest, and the overall *B* factor largest, indicating that the quality of the crystal grown in the solution was better than the others. Regarding the mosaicity, the average value (0.17) and deviation (0.06) were the lowest for the crystal grown in the PEG20,000 solution. The mosaicities of all of the frames are plotted in Figure 7. They showed that the anisotropy was smaller for the crystal grown from the PEG20,000 solution.

**Figure 7.** Mosaicities are plotted against the continuous frames: (**a**) 4834BK1 from PEG4000, (**b**) 4835M from PEG10,000 and (**c**) 4836M from PEG20,000.

In this study, lysozyme crystals with good diffractive quality were obtained in the viscous PEG solution using the counter-diffusion method. Their maximum resolutions were comparable to or better than the already registered data for a tetragonal lysozyme crystal in the PDB (for example, PDB ID: 4hp0, 5kxz, 6g8a), except for the crystal grown in space (1iee). Among the three MW PEGs, the best crystals were grown in PEG20,000. There are three possible reasons for this. (1) The density-driven flow around the growing crystals is reduced when the viscosity of the solution is higher. As a result, the highly ordered attachment of the protein molecules to the crystal surface suffers less from the flow [38,39]. (2) In a highly viscous solution, the reduction in the density-driven flow and the decrease of diffusive migration of the protein molecules enhances the formation of the protein depletion zone around the growing crystal [4]. As a result, the protein concentration on the surface of the growing crystal decreases, allowing for good crystals to be grown under lower supersaturation. (3) The excluded volume effect is stronger in high-MW-PEG solutions [40]. This generates a much stronger force to drive the protein molecules onto the crystal surface. The stronger macroscopic force on the surface of the crystal moves the protein molecules closer together inside the crystals and the molecular arrangement in the crystals becomes better. The mechanisms of growing good crystals in viscous high-MW-PEG solutions will be revealed by further investigations.

#### **4. Conclusions**

The diffusion mobility of solutes in viscous PEG solvents shows a considerable deviation, which is faster than the expected diffusion rate determined from the macroviscosity. This behavior was quantitatively described in this paper using an approximation model. The model enabled us to describe the mass transportation of molecules during the counterdiffusion crystallization processes. Figure 6a,b clarified that the self-searching mechanism for the optimal crystallization scenario would not work well when high-MW PEG is used as a precipitant in the counter-diffusion method.

However, crystallization in the counter-diffusion chamber containing a uniform concentration of PEG20,000 in advance works well. The greater extent of anomalous diffusion in a high-MW-PEG solution was observed with small-molecule solutes as salt precipitants, so that without a large delay in the crystallization period, the protein solution pre-mixed with viscous high-MW-PEG solutions could be crystalized using salt precipitants in the counter-diffusion crystallization method. Figure 6c shows that the self-searching mechanism worked well under this condition, and most of the protein-NaCl plane could be scanned after 16 days. We further confirmed the possibility of crystallization experiments with lysozyme solution pre-mixed with viscous high-MW PEG.

**Author Contributions:** Conceptualization, H.T. and A.T.; sample preparation and FCS experiments, R.U., K.T. and A.T.; crystallization experiment, A.T.; validation, A.T.; X-ray structural analysis, S.T. and H.T.; writing—original draft preparation, A.T.; writing—review and editing, A.T. and

H.T.; visualization, A.T. and H.T.; supervision, A.T. and H.T.; project administration, A.T.; funding acquisition, A.T. All authors have read and agreed to the published version of the manuscript.

**Funding:** This work was supported in part by the Development of Fundamental Technology for Protein Analyses and the Target Protein Research Programs from the Ministry of Education, Culture, Sports, Science and Technology of Japan.

**Institutional Review Board Statement:** Not applicable.

**Informed Consent Statement:** Not applicable.

**Data Availability Statement:** Not applicable.

**Acknowledgments:** We are grateful to Shigeyuki Yokoyama and Mikako Shirouzu (RIKEN) for their encouraging support.

**Conflicts of Interest:** The authors declare no conflict of interest.

#### **References**


### *Article* **The Effect of Controlled Mixing on ROY Polymorphism**

**Margot Van Nerom 1,\*, Pierre Gelin <sup>2</sup> , Mehrnaz Hashemiesfahan 2,3, Wim De Malsche <sup>2</sup> , James F. Lutsko <sup>4</sup> , Dominique Maes 1,\* and Quentin Galand <sup>1</sup>**


**Abstract:** We report the investigation of various experimental conditions and their influence on polymorphism of 5-methyl-2-[(2-nitrophenyl)amino]-3-thiophenecarbonitrile, commonly known as ROY. These conditions include an in-house-developed microfluidic chip with controlled mixing of parallel flows. We observed that different ROY concentrations and different solvent to antisolvent ratios naturally favored different polymorphs. Nonetheless, identical samples prepared with different mixing methods, such as rotation and magnetic stirring, consistently led to the formation of different polymorphs. A fourth parameter, namely the confinement of the sample, was also considered. Untangling all those parameters and their influences on polymorphism called for an experimental setup allowing all four to be controlled accurately. To that end, we developed a novel customized microfluidic setup allowing reproducible and controlled mixing conditions. Two parallel flows of antisolvent and ROY dissolved in solvent were infused into a transparent microchannel. Next, slow and progressive mixing could be obtained by molecular diffusion. Additionally, the microfluidic chip was equipped with a piezoceramic element, allowing the implementation of various mixing rates by acoustic mixing. With this device, we demonstrated the importance of parameters other than concentration on the polymorphism of ROY.

**Keywords:** ROY; polymorphism; microfluidics; acoustic mixing; diffusive mixing

#### **1. Introduction**

Polymorphism in crystallography was first defined by Eilhardt Mitscherlich in the 19th century and refers to the property of some chemical compositions to exist in multiple crystalline forms [1]. These different crystal structures are often a result of different molecular shapes caused by different torsion angles. Indeed, free rotation about single bonds within a molecule allows for several arrangements with potential energy minima and are therefore considered stable conformations. Consequently, the distinct molecular shapes result in different packing configurations. This particular form of polymorphism is termed conformational polymorphism [2]. However, different conformations exhibit different properties, such as stability, dissolvability, physiological activity, and/or bioavailability. This explains the importance of studies on crystal polymorphism in the pharmaceutical field. In drug development, the main objective is to achieve specific properties, and contamination with an undesired polymorph would be detrimental [3–6].

A model compound for studying polymorphism is 5-methyl-2-[(2-nitrophenyl)amino]- 3-thiophenecarbonitrile or ROY because of its ability to form at least 11 polymorphs, of which 6 are stable at room temperature. The latter are yellow prisms (Ys), orange needles (ONs), orange plates (OPs), red prisms (Rs), yellow needles (YNs), and orange-red plates (ORPs) and are formed by increasing the antisolvent-to-solvent ratio, which are water and

**Citation:** Van Nerom, M.; Gelin, P.; Hashemiesfahan, M.; De Malsche, W.; Lutsko, J.F.; Maes, D.; Galand, Q. The Effect of Controlled Mixing on ROY Polymorphism. *Crystals* **2022**, *12*, 577. https://doi.org/10.3390/ cryst12050577

Academic Editor: Abel Moreno

Received: 22 March 2022 Accepted: 14 April 2022 Published: 20 April 2022

**Publisher's Note:** MDPI stays neutral with regard to jurisdictional claims in published maps and institutional affiliations.

**Copyright:** © 2022 by the authors. Licensee MDPI, Basel, Switzerland. This article is an open access article distributed under the terms and conditions of the Creative Commons Attribution (CC BY) license (https:// creativecommons.org/licenses/by/ 4.0/).

acetone, respectively. The less stable and consequently less characterized forms are red plates (RPL), Y04, R05, YT04, and PO13. Red plates crystallize from vapor on succinic acid, Y04 and R05 can be obtained from a melt crystallization, and YT04 is a transformed form of Y04. PO13 is a supercooled melt form of YNs [7]. This abundance of polymorphs can be explained by the presence of three torsion angles and their different conformations. Moreover, the polymorphs are easily detected according to their color and shape, an aspect represented in the acronym Red Orange Yellow (ROY) [5–8].

Creating a controlled mixing regime that allows the selection of one particular polymorph with its specific properties is relevant for multiple industrial disciplines [3,9]. However, the crystallization of polymorphic compounds is a complex process and depends on many parameters, such as temperature, the solvent-to-antisolvent ratio, antisolvent properties, concentration, etc. As a first step toward a controlled and constant mixing environment, this work reports on the implementation of a microfluidic chip to achieve different controlled mixing regimes [10].

#### **2. Materials and Methods**

#### *2.1. Chemicals*

5-methyl-2-[(2-nitrophenyl)amino]-3-thiophenecarbonitrile (ROY) was purchased from TRC Canada® (Toronto, Canada) and consisted of a mixture of polymorphs with sizes ranging from about 1 to 200 µm (Figure 1). For stock solutions, synthesis-grade acetone from Acros Organics® (Geel, Belgium) with a purity of 99.9% and milliQ water were used. Stock solutions (starting solutions with ROY and antisolvent solutions without ROY) were prepared with a Sartorius CPA224S scale (Goettingen, Germany) with an accuracy of 0.1 mg. In the following, the concentrations of the samples are mentioned according to "acetone (Vol%)/Water (Vol%)/ROY (mg/mL)". The volumes of all prepared stock solutions were selected to be sufficiently large (>20 mL) so that the concentrations of all samples were reported with an accuracy of volume fractions below 0.01 Vol% and accuracy of ROY concentrations below 0.1 mg/mL.

**Figure 1.** ROY used in this study produced a mixture of polymorphs with sizes ranging from 1 to 200 µm.

#### *2.2. Crystallization in Bulk Experiments*

Antisolvent solution was added to and mixed with starting solution up to a total volume of 1.0 mL in 1.5 mL transparent vials at room temperature (RT). The respective volumes of the infused solutions were adjusted with KD Scientific® model 100 syringe pumps (Holliston, MA, USA) using Hamilton® Gastight® Luer Lock syringes (Allston, MA, USA). Connections between the syringes and the vials were realized with Idex connectors and 450 µm inner diameter (i.d.) and 670 µm outer diameter (o.d.) glass capillaries. The flow rate was set for each sample such that the total infusion time per sample was 10 min. During infusion, mixing was performed by rotating an 8 mm magnet at 100 RPM unless otherwise mentioned, or by shaking the vial with an Eppendorf thermomixer comfort at 500 RPM unless otherwise mentioned. Evaporation of the sample was prevented by sealing the glass vials (Figure 2). Each result documented in the following is based on at least three experimental runs unless otherwise mentioned.

**Figure 2.** Setup of crystallization in bulk experiments: (**a**) sample is mixed by a magnetic stirrer during 10 min infusion; (**b**) mixed by rotation during 10 min infusion.

#### *2.3. Crystzallization in Confinement*

To evaluate the possible influence of volume on the crystallization process of ROY, mixing was first performed as mentioned in the previous Section 2.2. Immediately after mixing and before crystallization occurred, the final solution was administered to capillaries with an i.d. of 100 µm, 200 µm, or 400 µm (Borosilicate glass round tubes from CM Scientific® (Silsden, UK)). The capillaries were then sealed with melted wax (Figure 3).

**Figure 3.** Confinement experiment in which sample is administered to capillaries before nucleation. Capillary is sealed afterwards with wax.

#### *2.4. Polymorph Identification*

Samples were analyzed with an Olympus® SZX16 microscope (Tokyo, Japan). The combination of the magnification of the microscope and the objective of the Olympus® SC100 color camera (Tokyo, Japan) resulted in a magnification of 23×. Different ROY polymorphs are so distinct in structure and color that optical identification was an easy and suitable method for this study. Moreover, experiments were repeated at least three times to fully characterize the observed sample.

#### *2.5. Microfluidic Setup*

All experiments were performed in microfluidic chips. The microchannel was etched in a silicon wafer, and the top and bottom of the channel were sealed with borosilicate glass. The internal section (width × depth) of the channels used in this study was 0.375 mm × 0.525 mm, and the channel length was 20 mm (Figure 4c,d). The channels were equipped with three inlets and outlets for infusion and evacuation of liquids and sealed by bonding 200 µm internal diameter Polymicro Technologies® (Phoenix, AZ, USA) glass capillaries with dual-cure epoxy sealant (Figure 4b). All connections were realized with IDEX connectors. Starting solution and antisolvent solution were simultaneously infused, and accurate flow rate controls were performed using Fluigent® pressure flow controllers (Lowell, MA, USA) (Figure 4a).

**Figure 4.** Microfluidic setup: (**a**) Simultaneous infusion of starting solution and antisolvent solution with accurate flow controls were performed with pressure flow controllers; (**b**) inlets and outlets were sealed with glass capillary tubes; (**c**) channel length of the chips was 20 mm; (**d**) internal section (width × depth) of the chips was 375 µm × 525 µm.

In-chip liquid mixing was performed by inducing acoustic streaming, for which a piezo-ceramic actuator (15 mm × 20 mm × 1 mm, APC international, Mackeyville, PA USA) with an eigenfrequency of around 2.0 MHz was placed at the back of the chip (Figure 5a). A Tektronix AFG1062 function generator was used to apply a sinusoidal voltage to the piezo element and an RF power amplifier (210 L, Electronics & Innovations, Rochester, NY USA) amplified the applied voltage with a maximal total output power of 10 W. To prevent the heating of the sample with the supplied acoustic energy, active temperature control was implemented: the sample was thermostabilized by a Peltier element driven by a PID controller with feedback from a thermistor incorporated in an aluminum block close to the working volume. A 120 W Peltier element was used, and a constant temperature water loop was applied at the backside of the thermoelectric element to improve the performance of the system (Figure 5b).

**Figure 5.** In-chip liquid mixing: (**a**) acoustic mixing of the liquid was obtained by vibrating the window of the chip with a piezo-ceramic actuator; (**b**) the cell system incorporated a PID temperature control system.

The overall setup was installed on an optical bench and imaged with a color HD camera through a 4× magnification objective. The microchannel was fixed on a translation stage, and the entire length of the channel could be placed in the field of view of the camera when searching for polymorphs. An overview of the experimental setup is pictured in Figure 6.

**Figure 6.** Overview of microfluidic experimental setup, including the PID temperature control system and the pressure flow controllers.

#### **3. Results**

#### *3.1. Supersaturation Protocol Influences ROY Polymorphism*

The most important parameters of crystallization of any compound are the concentration and the solvent-to-antisolvent ratio. ROY polymorphism has proven to be very sensitive to this ratio, and even the slightest change influences ROY polymorphism. To illustrate this, one batch shown here with a solvent volume percentage of 47.00 Vol% contained yellow needles (YNs) after mixing with a magnetic stirrer, while two batches with 46.00 Vol% and 45.00 Vol% solvent presented no YNs (Figure 7b). These minimal changes with different results prompted the implementation of repetitions between 3 and 20 for the following results.

To better control mixing, a setup was created in which flows were controlled so that different parameters could be tested. Within this setup, we observed that despite the same final position within the phase diagram of ROY, the supersaturation protocol to reach this position can alter ROY polymorphism. First, the initial concentrations of the starting solution and the antisolvent solution before mixing alter the final polymorphs observed after mixing. Indeed, mixing an antisolvent mixture of 100.00 Vol% antisolvent with a starting solution until a final volume percentage of 57.50 Vol% was reached resulted in orange plates (OPs) and orange-red plates (ORPs). Contrarily, mixing a solution of 70.00 Vol% antisolvent with a starting solution until the same volume percentage of 57.50 Vol% was reached resulted in yellow prisms (Ys) (Figure 7d,e).

Secondly, different mixing methods also alter the polymorph behavior independent of the final concentration and the solvent-to-antisolvent ratio (1.0 mg/mL ROY and 42.50 Vol% solvent, respectively). Mixing by a magnetic stirrer resulted in Ys, while mixing by rotation resulted in orange needles (ONs), yellow needles (YNs), and Ys (Figure 7e,f). Similarly, mixing by rotation until an end concentration of 5.0 mg/mL ROY was reached led to a different mix of polymorphs compared with samples mixed by a magnetic stirrer (Figure 7b). The effect of the mixing rate was also analyzed: the same set of polymorphs was consistently found, indicating that the mixing rate has no influence on polymorphism (Figures S1 and S2). In conclusion, the supersaturation protocol clearly affects ROY polymorphism. These observations led to the realization that a more controlled mixing method is required and thus motivated our implementation of several mixing methods in a microfluidic setup.

**Figure 7.** Supersaturation protocol influences ROY Polymorphism. (**a**) Schematic overview with possible trajectories of bulk mixing experiments performed within phase diagram; the supersaturation protocols represented here are mixing by rotation versus mixing by a magnetic stirrer. (**b**) Results after mixing with a magnetic stirrer (full line) or rotation (dotted line) with shown end concentrations and solvent to antisolvent ratios. (**c**) Schematic overview. Supersaturation protocols represented here are mixing by rotation versus mixing by a magnetic stirrer (blue versus red) and varying concentrations of the infused solutions (dotted versus full line). (**d**–**f**) Microscopic images of samples with the following concentrations and supersaturation protocol. (**d**) Mixed with a magnetic stirrer. Initial antisolvent solution: 100.00 Vol%. (**e**) Mixed with a magnetic stirrer. Initial antisolvent solution: 70.00 Vol%. (**f**) Mixed by rotation. Initial antisolvent solution: 100.00 Vol%. Unable to focus on all polymorphs simultaneously, the picture depicted here was carefully chosen with Ys in focus and ONs and YNs observed in the background.

#### *3.2. Confinement Has Minimal Effect on ROY Polymorphism*

Continuing in a microfluidic environment raises the question of whether confinement is also to be considered as a parameter that influences ROY polymorphism. To analyze this effect, samples were administered to capillaries with a diameter of 100 µm, 200 µm, or 400 µm before crystallization. The dimensions of the tested capillaries compare to the dimensions of the microfluidic chips. For all experimental conditions, the same polymorphs were found in all capillaries independently of the diameter. To mimic a more drastic change in volume, the comparison of the crystallization of identical samples between bulk and confinement was made. No major differences could be seen between bulk and confinement for the different samples tested, except for the following condition: 1.0 mg/mL ROY in 42.50 Vol% of solvent (Figure 8). In this sample, orange needles (ONs) with yellow needles (YNs) were observed when mixed in bulk, while orange plates (OPs) were observed when administered to a capillary. In order to identify the cause for these contrasting results, we set up a condition in bulk with the same solvent volume percentage and a lower ROY concentration of 0.8 mg/mL. In these conditions, OPs were formed, as they

were in the capillary with 1.0 mg/mL. This suggests that the effect of confinement on ROY polymorphism is minimal and can be explained by transient local depletion of ROY molecules that can nucleate within a confined space.

**Figure 8.** Confinement has minimal effect on ROY polymorphism. Obtained microscopic images of samples with concentrations and solvent to antisolvent ratios mentioned above the respective images. The upper part represents experiments performed in bulk. The bottom parts represent experiments in capillaries with an i.d. of 200 µm.

#### *3.3. In-Chip Mixing*

Two parallel flows of ROY and antisolvent solutions were infused simultaneously in the microchannel, and mixing was performed by inducing acoustic streaming. Acoustic streaming is known to be very sensitive to the experimental parameters. Indeed, acoustically soft water inside the channel surrounded by the acoustically hard silicon/glass chip forms an acoustic cavity. This implies that acoustic resonance occurs for certain specific frequencies. By tuning the applied frequency to one of the resonance frequencies so as to obtain a stationary wave along the width of the channel, the acoustic energy density inside the cavity is several orders of magnitude larger than it is at other frequencies [11], and the acoustic forces become strong enough to obtain efficient mixing [12]. The frequency of the applied signal was carefully adjusted to a resonance frequency with a precision of up to 0.001 MHz. Experimentally, this step was performed by infusing colored liquids and by monitoring the efficiency of the acoustic mixing as a function of the frequency. When an optimal frequency was reached, colored vortices were easily observed. Under these conditions, viscous attenuation of the acoustic wave in the liquid boundary layer resulted in four vortices in the direction of the acoustic propagation [11] (Figure 9). In addition, it is well known from the literature that the velocity of the liquid in the vortices increases approximately quadratically with the amplitude of the applied voltage [13]. Different mixing rates were obtained by varying the amplitude of the voltage applied to the piezoelectric element; voltage values of 0.12 V, 0.2 V, 0.5 V, and 1 V were used.

**Figure 9.** The frequency of the piezo element was tuned to the depth of the microfluidic channel to obtain stable vortices.

By applying sufficient acoustic power, complete mixing in the microfluidic channel occurred very rapidly. The amount of energy transferred was high, and this system caused the sample to heat up. In the absence of regulation, an increase of about 30 ◦C was observed within one minute. To counteract this, a PID temperature controller was used to stabilize the sample at room temperature. A typical stability of the measured temperature was ±0.1 ◦C (RMS value).

The two extreme cases tested are compared in Figure 10: in the first case, no voltage was applied, and the mixing occurred spontaneously by molecular diffusion. At time zero, when the flows were stopped, a clear interface between the liquids was observed. Mixing occurred gradually and was monitored qualitatively by the color variation. After about 60 s, no more color change was observed. In the second case, a 1 Vpp voltage was applied to the piezo element. Vortices instantaneously appeared, and complete mixing was obtained within about 2 s.

#### *3.4. Supersaturation Protocol Used within Microfluidic Channel Influences ROY Polymorphism*

A large number of experiments was performed in order to investigate the polymorphism in a microfluidic environment. Two main parameters were tested: on the one hand, different solvent and antisolvent concentration gradients were tested by varying the concentration of the infused solutions. For the first series of samples, depicted with dotted lines in Figure 11, pure antisolvent and ROY–acetone solutions were infused in the channel. For the other samples, depicted with solid lines in Figure 11, smaller solvent and antisolvent concentration gradients were created by infusing solutions with 70.00 and 55.00 Vol% antisolvent. Respective flowrates were adjusted so that after mixing, the final concentration in the microchannel corresponded to 42.50 Vol% acetone, 57.50 Vol% water, and 1.0 mg/mL ROY. These concentrations were identical to the concentrations of samples d, e, and f tested in bulk conditions in Figure 7.

Both the initial concentrations of the solutions and the mixing methods were found to play an important role in polymorphism. Indeed, large initial concentration gradients (Figure 11a) and diffusive mixing consistently resulted in yellow and orange needles (YNs and ONs), as shown in Figure 11b. The needles appeared in the central region of the channel, at the interface between the two liquids, and quickly developed in the entire volume.

However, with the same solutions, crystallization appeared to be very different when acoustic mixing was applied and mixtures of orange plates (OPs) and yellow prisms (Ys) were obtained for all tested samples, as shown in Figure 11c.

**Figure 11.** The supersaturation protocol influences ROY polymorphism in a microfluidic environment. Unable to focus on all polymorphs simultaneously, the pictures depicted here were carefully chosen with Ys in focus and the other polymorphs observed in the background. (**a**,**d**) Two supersaturation protocols were tested by varying the concentrations of the infused solutions. Influence of mixing rate was investigated by comparing diffusive mixing and acoustic mixing. (**b**) Larger initial concentration gradients and diffusive mixing consistently resulted in orange needles (ONs) and yellow needles (YNs). (**c**) Larger initial concentration gradients and acoustic mixing consistently resulted in orange plates (OPs) and yellow prisms (Ys). (**e**,**f**) Smaller concentration gradients and diffusive mixing produced different results: in 7 samples out of 10, orange plates (OPs) and yellow prims (Ys) were obtained. In 3 samples, yellow needles (YNs) were also observed. (**g**) Smaller initial concentration gradients and acoustic mixing consistently resulted in orange plates (OPs) and yellow prisms (Ys).

With smaller concentration gradients (Figure 11d), diffusive mixing resulted in different polymorphs. For most samples (7 out of 10), a mixture of OPs and Ys was obtained (Figure 11f), while for three samples, large YNs were also observed (Figure 11e). Orange plates were mostly observed along the walls of the channel, but crystals appeared everywhere. Finally, for those solutions and in all tested samples, acoustic mixing resulted in OPs and Ys, as shown in Figure 11g.

In these samples, crystallization began as soon as mixing started along the walls of the microchannels, and crystals were observed in the entire volume in a few seconds. In most experimental runs, OPs grew faster, and Ys were detected 30 s to one minute after the mixing. Identical polymorphs were obtained for all tested intensities of acoustic mixing (Figure S3).

#### **4. Discussion**

From the overall experimental dataset obtained in this study, we made some general key observations about ROY polymorphism. ROY polymorphism is very sensitive to the solvent-to-antisolvent ratio, as shown in the results obtained in Figure 7. Moreover, it appears that mixing is a key parameter. A homogeneous sample after mixing by rotation is reached after several seconds. On the timescale of seconds, different local concentrations are possible in the vial during the entire infusion process, and strong local concentration gradients exist in the mixture. This allows the system to sample several regions in the free energy diagram and nuclei of various polymorphs to form and reach critical sizes. This

results in very complex final mixtures of many polymorphs, as shown in the dotted lined samples in Figure 7b. This statement is supported by the results displayed in Figure 7d,e, in which different start concentrations were used to mix to the final end concentration. Indeed, lowering the initial concentration gradients results in smaller local concentration gradients, resulting in a smaller variety of polymorphs in the end. Contrarily, mixing with a magnetic stirrer allows faster mixing and keeps the concentration homogeneous in the entire vial. This mixing method results in less drastic local concentration gradients, prohibiting the system from reaching certain regions in the free energy diagram. This gives rise to a smaller variety of polymorphs compared with the samples mixed by rotation.

Similar observations emerged when analyzing the results of microfluidic experiments. The results in Figure 8 show that the influence of the confinement of the sample is low. However, for identical final conditions, two mixing methods within the microfluidic channels led to other polymorphs, as shown in Figure 11. Diffusive mixing is slower, and the concentration gradients in the liquid result in less reproducible experiments allowing the formation of different polymorphs, while acoustic mixing always leads to orange plates and yellow prisms. Acoustic mixing is very efficient, and the system reaches its final mean concentration in the entire microfluidic channel within 1 to 2 s. These results indicate that in the final concentration condition, that is, 47.50 Vol% acetone, 52.50 Vol% water, and 1.0 mg/mL ROY; orange plates and yellow prisms are the most stable polymorphs.

It also appears that needles can be obtained with rotative mixing in bulk or with diffusive mixing in a microchannel, but over the very large number of experiments that we performed, needles were never obtained with magnetic stirring in vials nor with acoustic mixing in microchannels. These observations suggest that needle nucleation occurs on longer timescales. On the basis of our experiments, however, it is not possible to draw a solid conclusion, and the formation of needles could also be hindered by excessive shear and flow resulting from the magnetic stirring and acoustic mixing. Moreover, a separate study on nucleation may bring interesting insights into the driver of these differences observed.

Throughout all our experiments, orange-red plates were observed in very few samples. This polymorph is not stable in the concentration domain that we tested. Finally, it appears that the microfluidic environment and the use of fast mixing methods represent an interesting opportunity for the selective production of certain polymorphs. Under the conditions that we tested, the preparation of the samples in microfluidic reactors and using acoustic mixing consistently produced orange plates and yellow prisms.

#### **5. Conclusions**

This work reports on ROY crystallization in various conditions. Not only do the final end concentration of the solvent, antisolvent, and ROY affect polymorphism, but the selected protocol used to reach supersaturation also affects the polymorphism behavior. Here, we reported on two different sets of protocols resulting in different polymorphs: first by changing the initial mixing solutions, and second, by altering the mixing method. In bulk, mixing was performed either by rotation or a magnetic stirrer. In a microfluidic setup, we created controlled mixing environments and tested the two sets of protocols: changing the initial concentrations and the mixing method. In the microfluidic channel, this was achieved by either diffusive mixing or mixing acoustic streaming. To analyze mixing in a microfluidic setup, the influence of confinement on polymorphism was considered. The only effect of confinement that we discovered was the transient local depletion of ROY molecules.

Both in bulk and in a microfluidic environment, changing the initial stock solutions resulted in different polymorphs, as did the applied mixing method. We believe that mixing by rotation in bulk and mixing by diffusion in a microfluidic channel allows the system to reach minima in the free energy diagram because of slow mixing and large local gradients, which are not reached (or do not last long enough) when fast mixing is performed by a magnetic stirrer in bulk or acoustic streaming in a microfluidic channel. This results in more variable and less stable polymorphs in the conditions tested when slow mixing occurs. On

the basis of these findings, we conclude that orange plates and yellow prisms are the most stable polymorphs in the system that we tested. Indeed, needles only appear under gentle mixing conditions, such as diffusion and mixing by rotation, while orange-red plates were sporadically obtained.

In the microfluidic setup, in which mixing can be accurately controlled, orange plates and yellow prisms can be efficiently and reproducibly selected. Overall, these findings show that the start conditions, supersaturation protocol, and flow are critical parameters in ROY polymorphism that have to be considered in the design of systems for polymorph selectivity.

**Supplementary Materials:** The following supporting information can be downloaded at: https: //www.mdpi.com/article/10.3390/cryst12050577/s1, Figure S1: mixing rate in bulk experiments with rotative mixing does not influence polymorphism; Figure S2: mixing rate in bulk experiments with magnetic stirrer does not influence polymorphism; Figure S3: intensity of acoustic mixing does not influence polymorphism in microfluidics.

**Author Contributions:** Conceptualization, D.M., W.D.M., Q.G. and M.V.N.; methodology, D.M., Q.G., M.H. and M.V.N.; experiment design, Q.G., M.V.N., M.H. and P.G.; resources, M.H.; writing—original draft preparation, M.V.N. and Q.G.; writing—review and editing, M.V.N., Q.G., D.M. and J.F.L.; supervision, D.M. and W.D.M.; project administration, D.M.; funding acquisition, D.M. All authors have read and agreed to the published version of the manuscript.

**Funding:** This research was funded by the European Space Agency under Prodex Contract No. ESA AO-2004-070 and by a Strategic Research Program on Microfluidics (SRP51) at Vrije Universiteit Brussel.

**Institutional Review Board Statement:** Not applicable.

**Informed Consent Statement:** Not applicable.

**Data Availability Statement:** Not applicable.

**Conflicts of Interest:** The authors declare no conflict of interest. The funders had no role in the design of the study; in the collection, analyses, or interpretation of data; in the writing of the manuscript; or in the decision to publish the results.

#### **References**


**Benjamin Radel \*, Marco Gleiß and Hermann Nirschl**

Institute of Mechanical Process Engineering and Mechanics, Karlsruhe Institute of Technology, 8 Strasse am Forum, 76137 Karlsruhe, Germany; marco.gleiss@kit.edu (M.G.); hermann.nirschl@kit.edu (H.N.) **\*** Correspondence: benjamin.radel@kit.edu

**Abstract:** Combined normal and shear stress on particles occurs in many devices for solid–liquid separation. Protein crystals are much more fragile compared to conventional crystals because of their high water content. Therefore, unwanted crystal breakage is to be expected in the processing of such materials. The influence of pressure and shearing has been investigated individually in the past. To analyze the influence of combined shear and normal stress on protein crystals, a modified shear cell for a ring shear tester is used. This device allows one to accurately vary the normal and shear stress on moist crystals in a saturated particle bed. Analyzing the protein crystals in a moist state is important because the mechanical properties change significantly after drying. The results show a big influence of the applied normal stress on crystal breakage while shearing. Higher normal loading leads to a much bigger comminution. The shear velocity, however, has a comparatively negligible influence.

**Keywords:** protein crystals; breakage; shear stress

**Citation:** Radel, B.; Gleiß, M.; Nirschl, H. Crystal Breakage Due to Combined Normal and Shear Loading. *Crystals* **2022**, *12*, 644. https://doi.org/10.3390/ cryst12050644

Academic Editor: Abel Moreno

Received: 31 March 2022 Accepted: 27 April 2022 Published: 30 April 2022

**Publisher's Note:** MDPI stays neutral with regard to jurisdictional claims in published maps and institutional affiliations.

**Copyright:** © 2022 by the authors. Licensee MDPI, Basel, Switzerland. This article is an open access article distributed under the terms and conditions of the Creative Commons Attribution (CC BY) license (https:// creativecommons.org/licenses/by/ 4.0/).

#### **1. Introduction**

In recent years, alternative methods to purify and formulate proteins have been investigated. One example is preparative protein crystallization [1], which provides the opportunity to replace a costly chromatography step with selective protein crystallization. Crystalline proteins offer several advantages. By influencing the crystal shape and size, product characteristics like handling, shelf life, and drug release properties can be adjusted [2,3]. One well-known example is the use of crystalline insulin to achieve a depot effect and constant bioavailability [4]. Conventionally, protein crystallization on a larger scale uses the displacement method. Added anti-solvents (typically salts) reduce the target protein's solubility and, thus, create the required supersaturation. Inherently, this method has some major disadvantages. First, the anti-solvent is typically added as a solution, which dilutes and, thus reduces the final crystal concentration in the suspension. Second, the final solution has high ion strength and, finally, high local supersaturation at the anti-solvent inlet leads to a broad crystal size distribution, lower reproducibility, and increased danger of amorphous precipitation [5]. Evaporative crystallization tackles these disadvantages, but is not suitable for proteins because of the required high temperatures. Groß and Kind [5] introduced low temperature water evaporation crystallization for proteins and demonstrated the applicability for the model system lysozyme from hen-egg white. Small scale filtration experiments of such crystallizate by Radel et al. [6,7] allowed to determine the process functions of this particle system. Using these process functions, the low temperature water evaporation crystallization method was adapted and implemented by Dobler et al. [8] for an integrated, quasi-continuous apparatus concept. Barros Groß and Kind [9] further investigated how seeding affects the crystal size distribution and demonstrated high reproducibility and control over crystal sizes and the width of the size distribution.

One major difference between protein and conventional crystals is the high water content in the protein crystal, which can account for up to 80% of the crystal mass. This makes such crystals soft and sensitive to mechanical stress [10]. Cornehl et al. [11] investigated crystal breakage of aggregated and needle-shaped lysozyme crystals under compressive stress. In process engineering, compressive stress occurs, for example, in press or cake filtration. Other apparatuses for filtration, such as the cross-flow filter (BoCross Dynamic, BOKELA GmbH, Karlsruhe, Germany) are characterized by high shear stress due to integrated agitators. In these cases, Cornehl et al. [12] also observed crystal comminution. In crystallization itself, comminution is partly a desirable effect. For example, collision of the crystals with each other as well as with internals, walls, and stirrers creates new crystallization nuclei that subsequently grow into larger crystals. For the use of crystallization as a formulation step for downstream processing of proteins [3,13], comminution after crystallization is usually not desired, since a change in particle size distribution (PSD) also changes the product properties.

Sediments can be compacted by both normal and shear stress. If both types of stress are combined, the achievable compaction is more pronounced, as used, for example, by Illies et al. [14] for the dewatering of filter cakes. Höfgen et al. [15] used an apparatus with high pressure dewatering rolls for filtration and dewatering. With this system, the shear and normal forces are individually adjustable and adaptable to the application at hand. Hammerich et al. [16] modified a shear cell for the Schulze ring shear tester to allow the measurement of fluid saturated particulate networks. This provides the opportunity to study the rheology and flow behavior of saturated sediments. The higher achievable compaction with combined shear and compression indicates bigger mechanical stress. Hence, a more pronounced comminution for protein crystals is to be expected in this case. The setup with the modified shear cell can, therefore, be used to observe comminution at a defined normal stress and shear velocity.

#### **2. Theory**

Particle breakage occurs when the material strength is exceeded. In addition to a one-time high mechanical load, fracture due to several lower load cycles is also conceivable. In this case, the loading history of the particle must be taken into account. For comminution processes such as grinding, the occurrence of particle or aggregate breakage is the basic requirement. In most solid–liquid separation applications, particle breakage is undesirable and usually leads to a deterioration of the process result. Basically, a distinction must be made between the stress in a particulate network and on the individual particle. In the following, the special focus is on compressive and shear stress.

In a particulate network under compressive stress, the imposed mechanical load is degradable by rearrangement processes, deformation, and fracture of the particles. For an axial, one-dimensional load, such as that applied by a piston, the effects mentioned above result in the highest load being applied in the immediate vicinity of the piston. The load decreases with increasing distance from the piston in the sediment. If the sediment is compressible, the highest compaction is to be expected in the vicinity of the piston [17]. In the upper particle layers, absorption of the input mechanical energy takes place due to elastic and plastic deformation and particle breakage, such that the mechanical energy transmitted via contact points decreases. Thus, the lowest stresses are present at the greatest distance from the piston center and at the bottom [18].

Shear stresses occur in process engineering, for example, on agitators, valves, and in pumps. Fluid-induced mechanical stress occurs in a suspension due to the existing flow. In addition to the fluid-induced stress, collisions between particles and particles with agitators or walls take place. These collisions lead to particle abrasion and, at high loads, to particle breakage. In crystallization, such abrasion is partly desirable, as it creates new crystallization nuclei for secondary nucleation.

Shear loading is also possible for saturated sediments. This applies to the transport of flowable sediments in equipment, such as the tubular centrifuge or decanter centrifuge, as well as to the targeted use of combined normal and shear loading for mechanical dewatering. When normal and shear forces are superimposed, significantly higher densities and, thus,

lower residual moisture content can be achieved [14,15]. In such scenarios, the particles rub against each other and, in the case of mechanically labile particle systems, abrasion or fracture occurs.

#### *Comminution*

The PSD, particle shape, and sediment structure have a big influence on the energy absorption and stress distribution in the particle bed. Furthermore, the overall stress is affected by the normal and shear stress as well as the shear velocity. As stress increases, comminution occurs in addition to compaction. Comminution leads to new fracture surfaces and rearrangement processes. The smaller the particles, the higher the energies required to cause comminution.

Population balances can be used to model the comminution of different particle size classes. The probability of breakage, *P*, during impact loading is given for many materials by the exponential function

$$P = 1 - \exp\left(-f\text{x}k(E\_M - E\_{M,\text{min}})\right) \tag{1}$$

with *f* as the material parameter, *x* as the particle diameter, *k* as the number of load cycles, *E<sup>M</sup>* as mass-related stress energy, and *EM*,min as the threshold value of the mass-related stress energy [18]. The breakage fraction *P*¯ is identical to the breakage probability for a single grain. In a collective, the breakage fraction can be described by normalization with the master curve

$$\frac{\tilde{P}}{P\_{\infty}} = 1 - \exp\left(-\left(\frac{E\_M}{E\_{M,c}}\right)^{\beta}\right). \tag{2}$$

The quantity *P*¯<sup>∞</sup> is the upper limit of *P*¯ and *β* is a curve parameter. *EM*,*<sup>c</sup>* is a characteristic value of *E<sup>M</sup>* for which different approaches, taking into account the particle size, exist [19].

A fracture force can be determined on a single grain via compression or indentation tests. Nanoindentation is a suitable method for small particles that are to be examined moist, as in the case of protein crystals. Depending on the particle shape, different loads are conceivable. For an elongated particle, for example, a three-point bending test is possible. To cause particle fracture in a particulate network, a higher force than the fracture force of the single particle is necessary. This is due to the absorption of mechanical stress by rearrangement and deformation processes. Particle breakage is divided into the following phases:


Cracking starts as soon as the material strength is exceeded. Therefore, cracks preferentially form at locations that are already under stress or where defects are present. Defects can be, for example, lattice defects in a crystal. When the energy is low, a dormant crack is formed in this way. If no further stress takes place, the crack does not propagate and no fracture occurs. Only at a sufficiently high energy does crack initiation and subsequent crack propagation occur. The energy required for crack propagation at the cracking front must be continuously replenished.

Therefore, two conditions apply to particle breakage. On the one hand, there is the force condition, i.e., overcoming the binding forces for the formation of a crack, and on the other hand, there is the energy condition. The latter states that energy consumed at the crack front must be continuously supplied [20].

#### **3. Materials and Methods**

#### *3.1. Crystallization*

Isometric lysozyme crystals are produced with displacement crystallization. Two stock solutions are prepared. Solution one is an 25 mmol L−1 acetic acid buffer at pH 4.5. Solution two has the same composition but also contains 80 g L−1 NaCl. After dissolving 100 g L−1 lysozyme (Granulated lysozyme, OVOBEST Eiprodukte GmbH & Co. KG, Neuenkirchen-Vörden, Germany) from hen egg-white in 125 mL of stock solution one, 125 mL of stock solution two is added at a rate of 1 mL min−1 with a membrane pump. Afterwards, the crystals grow in the aging phase until the supersaturation is reduced to zero. During the whole crystallization process, which takes about 16 h, the solution is stirred with a blade stirrer at 350 min−1 .

#### *3.2. Modified Shear Cell and Ring Shear Tester*

For the combined normal and shear loading, sediments of isometric lysozyme crystals are used as a model system for protein crystals. The mechanical properties of dried and moist crystals differ greatly. Therefore, for a realistic assessment of the material behavior, it is important to measure the crystals as close as possible to the actual process conditions. To load the crystals in a saturated sediment, the ring shear tester (RST-01.pc, Dr. Dietmar Schulze) is used. This device is common in the field of bulk mechanics to characterize the flow properties of dry powders. The bulk material is placed in a shear cell. The shear cell lid with attached drivers is placed on the bulk material and hooked onto a weight. The normal stress of the system is determined by the applied weight. For shearing, the shear cell rotates, but the shear cell lid is connected and fixed to bending beams via tie rods. The shear stress can, therefore, be measured at the bending beams.

Modifications to the shear cell are necessary for the measurement of moist sediments. A detailed explanation of the modified shear cell can be found in Hammerich et al. [16]. Therefore, the modifications are only roughly outlined in the following. When moist sediments are loaded, compaction results in the displacement of fluid, which must escape from the system. In contrast, when the sediment is stretched back, it must be ensured that no desaturation occurs.

The structural implementation of the ring-shaped modified shear cell is shown as a sectional view in Figure 1. The shear cell lid consists of two parts, which are sealed against each other by Teflon rings. To allow liquid to escape, but at the same time retain the particles, filter media with a support structure are installed on the top and bottom of the shear cell, shown in green in Figure 1. The drivers and the bottom plate of the shear cell are visualized in gray. On the bottom of the shear cell, below the filter medium, a drainage channel is located, which allows the displaced fluid to escape. The L-shaped profile ensures that the drainage channel is always filled with liquid. Thus, the sediment does not desaturate in the event of back expansion. For the same purpose, there is a riser tube for the displaced liquid on the shear cell lid.

Due to the additional Teflon seals and the slightly different design, the measured values for normal and shear stress cannot be determined directly. A correction for the friction contributions of the seals is required. For this purpose, Hammerich et al. [16] developed a variant of the shear cell with strain gauges for load measurement and determined correction functions for the normal and shear stress. These functions can be found in [16]. The correction functions depend on the seal set used and, therefore, must be recalibrated when replacing the seals.

An experiment with the ring shear tester and the modified shear cell is divided into the following phases:


**Figure 1.** Sectional view of the modified shear cell; adapted from [16].

Sample preparation consists of preparing isometric lysozyme crystals via displacement crystallization. To provide a sufficient amount of crystallizate, and to compensate for minor variations between the different batches, three individual crystallization batches have to be prepared and mixed. The crystals settle overnight. The following day, the supernatant can be decanted. This increases the concentration of crystallizate when the sediment is resuspended and the higher solids volume fraction prevents segregation due to zone sedimentation.

When assembling the shear cell, the filter medium (Trakedge 0.2 µm, Sabeu GmbH & Co. KG, Northeim, Germany) the support structures, and the base plate have to be installed in the bottom part of the shear cell. The shear cell bottom is placed into a customized centrifuge insert for the beaker centrifuge (ZK 630, Berthold Hermle AG, Gosheim, Germany). This insert allows the sediment to be built up by centrifugation directly in the shear cell. Suspension is added to the beaker, so that the processing chamber of the shear cell bottom is slightly overfilled. The sample is centrifuged at a speed of 1500 min−1 for 10 min. This corresponds to a maximum relative centrifugal acceleration of *C* = 500. Due to the comparably low centrifugal acceleration, the normal loading of the particulate network is low and a change in PSD due to centrifugation can be avoided. After removal from the centrifuge, the sediment has been built up in the shear cell processing chamber. The supernatant is discarded carefully. The L-shaped drainage channels are filled with supernatant. To ensure a defined sediment level for the different tests, excess sediment is carefully removed using a scraper with a fixed length.

The same filter membrane, support structures, and drivers must also be mounted on the shear cell lid. The shear cell lid consists of two parts and, therefore, requires sealing with Teflon O-rings. The assembled lid is centered and placed on the shear cell with the sediment. By softly pressing on the lid, it is fixed in the correct orientation.

The shear cell preparation takes place next. The assembled shear cell is inserted into the ring shear tester device, the tie rods are installed, and the counterweight is attached to the lid. To finally position the shear cell lid in a defined way, a normal load of 40,180 Pa is applied. This causes the lid to slide into position and displace some fluid from the sediment. The riser tubes are then filled with supernatant. Afterwards, the normal load is removed, and the shear cell is now ready for the actual shear test.

The shear test consists of applying the defined normal stress and then shearing the specimen by rotating the bottom part of the shear cell. The adjustable parameters are the shear velocity, the normal stress, and the shear duration. For the characterization of powders, the default shear velocity is 1.5 mm min−1. The velocity is varied in the range of 0.48–4.5 mm min−1 for the combined shear and normal loading of the moist protein crystals. In the experiments, the shearing time is varied in such a way that the shearing path remains constant. A schematic plot of the shear stress versus shear duration is shown in Figure 2. After removal of the sheared sediment from the shear cell, samples for measuring the PSD

with laser scattering are taken and resuspended in supernatant. Additionaly, samples are cut out of the sediment with a polyimide tube for characterization with micro computed tomography (µCT). To assess the influence of centrifugation and the compaction at shear cell preparation on the PSD, samples of the initial solution and after centrifugation and compaction without shearing in the ring shear apparatus are also analyzed.

**Figure 2.** Schematic representation of the shear stress curve.

#### *3.3. Micro Computed Tomography Analytic*

The working principle and sample preparation as well as post processing are explained in detail in Dobler et al. [8]. Thus, the preparation of the samples is only briefly outlined in the following. The polyimide tubes holding the samples of the compressed and sheared sediment are deep-frozen with liquid nitrogen. Afterwards, the pore fluid is removed with lyophilization. The removal of pore fluid ensures better contrast between particles and the surrounding air. For the measurement in the µCT (Zeiss Xradia Versa 520, Carl Zeiss Microscopy GmbH, Oberkochen, Germany), the sample tube has to be glued onto a dress pin and mounted on a sample holder. At an energy of 50 W, 2201 projections (X-ray images) are taken from the sample at different rotation angles. These projections are reconstructed to a 3D 16 grayscale image stack. The resulting voxel size is about 400 nm.

#### **4. Results and Discussion**

The shear stress applied in the shear cell has a gradient and, thus, varies in the sample. Directly at the shear cell lid, between the drivers protruding into the sample, the shear stress is zero. Below the drivers, the highest load occurs, which then decreases towards the bottom of the shear cell. This also results in differences with regard to the comminution that occurs. No comminution takes place directly between the drivers of the lid.

The resulting shear zones are also visible in µCT images of the loaded sediment. Figure 3 shows an image of the top of the sediment after shearing with a normal stress of 88,396 Pa. No comminution, but only compression, takes place in the uppermost layer of the sediment (green). This is the area between the approximately 2 mm long drivers, which protrude into the sample. When disassembling the shear cell, a small amount of the sediment sticks to the shear cell lid. Thus, the height of the undamaged crystal layer is approximately 600 µm. The impression of these drivers in the sediment are also visible on the top of the µCT images in Figure 4. The particles in the area between the grooves of the top layer in Figure 4 experience no shear force. In Figure 3, the sharply delineated, red-colored layer directly below has the highest shear force. The resulting comminution of

the crystals is so strong, that the crystal structure can no longer be visualized with the µCT. The particle sizes and the porosity is too low for the CT scan resolution. This is shown by the fact that practically no crystal edges can be seen.

**Figure 3.** Shear zones in the sediment: **Top** (green), undamaged crystals in the area between the drivers; **bottom** (red), crushed and no longer resolvable crystals.

**Figure 4.** Computer tomography image of a bulk (gray) with tracer particles (green) in the initial unloaded (**left**) and sheared (**right**) state.

Comminution at the bottom of the shear cell and at the very top between the drivers is much less pronounced than at the top, directly below the drivers, of the saturated sediment. This fact is also confirmed by µCT scans of a model sediment in the initial unloaded and sheared condition in a miniaturized shear cell. The visualization of these scans is shown

in Figure 4. The particles colored in green are a model system with lower density, which therefore have a different gray value in the tomography. Thus, the initial and final position of these particles can be marked. The grooves seen on the top of the sediment are from the drivers. To better visualize the green particles, the other gray material in the front area is shown faded. In the initial state on the left side of Figure 4, the green tracer particles are present in a vertical stripe in the sediment. After a shear path of about 45°, the distribution of the tracer particles shown on the right is obtained. A large part of the tracer particles remains unsheared and stays together in a strip. In the upper region near the lid, however, there is a clear distribution of the tracer material along the shear path. In this area, the particles experience high shear stress, which, in the case of mechanically sensitive materials, also causes comminution. The distribution of the particles also explains why not all large particles are crushed. In the lower area of the sediment and at the top between the drivers, basically no comminution takes place.

Figure 5 shows the PSD of the unloaded suspension and of two layers of the sediment after shearing. The comminution of the sheared samples compared to the unloaded initial suspension can be seen. Furthermore, it can be seen from Figure 5, that the PSD of the top layer of the bulk has a higher fines content than the lower layer and the total sediment. The sheared total curve in Figure 5 is a laser diffraction measurement of a resuspendend sample of the sediment with complete sediment height. This curve, therefore, represents all layers in the sample. The comminution is particularly pronounced for the *x*<sup>10</sup> diameter, which is only between 10–15% of the unloaded suspension. The *x*<sup>90</sup> diameter decreases much less and is in the range between 75–85%. This can be explained by the aforementioned unequal shear stress in the sediment. Large particles remain intact at the bottom of the bulk and at the same time an increase in very small particles is observed. Thus, the PSD broadens.

**Figure 5.** Particle size distributions in the different bulk layers after shearing.

Various variables influencing comminution are conceivable. In particular, the applied normal stress and the shear velocity probably have an effect on crystal breakage. These two influencing variables are, therefore, considered in more detail in this study. The covered shear path for all experiments remains a constant 15 mm in order to exclude the influence of different shear paths. For the evaluation of the changes of the PSD, the comparison between the respective unloaded suspension and a resuspended sample of the total sediment after loading is used. The PSD is obtained with laser diffraction.

#### *4.1. Influence of Centrifugation and Compaction*

It is conceivable that centrifugation and compaction for sample preparation in the ring shear tester also causes particle size reduction. To assess this influence, the PSDs of unloaded and centrifuged and compacted samples are shown in Figure 6. Compacted in this context means an applied normal load of the sediment in the ring shear tester with *σ*<sup>n</sup> = 40,180 Pa, without shear loading and after centrifugation.

**Figure 6.** Comparison of particle size distributions for the unloaded suspension and centrifuged and compacted sediments.

The resulting PSDs can be seen in Figure 6. The crystals in Figure 6 are centrifuged at a speed of 1500 min−1 and subsequently compacted in the ring shear tester. Afterwards, the sediment has been resuspended in supernatant for the determination of the PSD with laser diffraction. They show a slight shift towards smaller particle sizes. However, the whole PSD is shifted and no increase in fines is observed. The shift is so weakly pronounced that comminution is negligible compared to the sheared samples. This slight change in PSD can also be explained by minor temperature variation during measurement and is not necessarily due to centrifugation. Hence, the preparation of the shear cell and the subsequent sampling, does not affect the PSD significantly. In summary, therefore, neither the load in the centrifuge nor the additional normal load in the shear cell causes significant particle size reduction or crystal breakage.

#### *4.2. Influence of Normal Stress*

Shear tests at a constant shear velocity of 1.5 mm min−1 and varying normal stresses are used to evaluate the influence of the applied normal stress on comminution. The transformed particle size density distribution *q* ∗ 3 for varying normal stress is shown to visualize the comminution. Since this is a biological system that has natural fluctuations during crystallization, the initial distribution is also displayed for each normal stress.

The resulting PSDs are shown in Figure 7. Already at a low normal load of 20,572 Pa a size reduction, as well as a shift of the particle sizes, can be observed and the distribution width increases by about 40%. The *x*<sup>90</sup> and *x*<sup>50</sup> diameters decrease significantly at this low load to 78 and 68% compared to the unloaded particle sizes. The most significant decrease is observed for the *x*<sup>10</sup> diameter. This diameter is reduced to 30% after loading.

**Figure 7.** Influence of different normal stresses on the particle size density distribution.

At a higher normal stress of 41,787 Pa, the same effects occur more strongly. The distribution width increases by 62% compared to the unloaded specimen. The characteristic diameters *x*<sup>90</sup> and *x*<sup>50</sup> are now only 74 and 58% of the initial sizes, respectively. The *x*<sup>10</sup> diameter also decreases significantly to 19%. The density distribution shows an increase in fines and a reduction of large particles. Compared to the normal stress of 20,572 Pa, the comminution at 41,787 Pa is more pronounced and the PSD is significantly wider. In the density distribution, this fact is shown by the reduction and broadening of the peak.

A further increase of the normal stress to 64,288 Pa leads to an obvious change in the particle size density distribution. The peak, still clearly visible in the unloaded sample, has almost completely disappeared. A strong increase of very small particles and an intense broadening of the distribution can be observed. The amount of larger particles decreases somewhat, but is still present in the shear-loaded sample. The distribution width is 147% higher compared to the unloaded sample. The *x*<sup>90</sup> diameter decreases to 63%. The median diameter also decreases significantly and is 32% of the initial diameter. This is a much greater decrease compared to the lower normal loadings. The *x*<sup>10</sup> diameter amounts to 15% of the initial unloaded diameter.

An increase of the normal stress beyond this to 88,396 Pa does not cause any further significant comminution of the particles. The distribution width is broadened with an increase of 136%. This is similar to the previous normal stress. The *x*<sup>90</sup> diameter is 68% of the initial diameter. The median and *x*<sup>10</sup> diameters are reduced to 38 and 12%, respectively. These values are also similar to the previous normal stress level. There is a clear increase in small particles.

The comparison of the different normal stresses shows a clear influence of the applied stress on the occurring comminution. Low normal stresses result in less comminution. However, even small normal stresses cause a widening of the PSD and a significant decrease of the *x*<sup>10</sup> diameter. With an increase of the normal stress, the comminution intensifies and the distribution width increases significantly. The median diameter is strongly reduced. Since the shear load in the sediment has a gradient, large particles are always retained, which is reflected in a less pronounced reduction of the *x*<sup>90</sup> diameter. Above a normal stress of 64,288 Pa, an increase in normal stress does not result in additional comminution. The decrease in characteristic diameter, increase in distribution width, and comminution are similar. This indicates an upper limit above which the normal stress loses its influence. The comminution is probably not only due to particle breakage but also to abrasion. Abrasion is caused by the friction of the crystals against each other, which, for example, abrades the corners of the crystals.

#### *4.3. Influence of Shear Velocity*

In addition to the normal stress, an influence of the shear velocity on the comminution is also conceivable. In order to investigate this influence, sediments of isometric lysozyme crystals are loaded with the three shear velocities 0.48; 1.5 and 4.5 mm min−1 at a constant normal stress of 64,288 Pa and a constant shear path. The particle size density distributions of the unloaded and sheared specimens with different shear velocities are shown in Figure 8.

Even at the lowest shear velocity in Figure 8, top left, there is a strong comminution of the particle collective. The distribution width increases by 138% and the *x*<sup>90</sup> and median diameters decrease to 75 and 41%, respectively. The *x*<sup>10</sup> diameter reduces significantly to 15% compared to the unloaded suspension. Increasing the shear velocity to 1.5 mm min−1 shows higher comminution for the *x*<sup>90</sup> and *x*<sup>50</sup> diameters and an increase in the distribution width. The respective values are analogous to Figure 7 at the corresponding normal load of 64,288 Pa. A further increase of the shear velocity to 4.5 mm s−1 leads to a very similar particle size density distribution of the sheared crystals compared to the previous shear velocity. The reduction of the characteristic median and *x*<sup>90</sup> diameter is more pronounced by a few percentage points.

The lower right plot in Figure 8 shows the particle size density distribution of a sediment after shearing at a speed of 4.5 mm min−1 and a normal stress of *σ<sup>n</sup>* = 20,572 Pa. Here, it can be seen that the normal stress has a stronger influence on the comminution than the shear velocity. At a high shear velocity and low normal stress, the comminution is less. The *x*<sup>90</sup> diameter is still 80% of the initial particle size. The median diameter decreases to 69% and the *x*<sup>10</sup> diameter reduces significantly to 15%. The reduction in the median and *x*<sup>90</sup> diameters thus corresponds almost exactly to the values at a normal load of 20,572 Pa and a shear velocity of 1.5 mm min−1 from Figure 7. However, the reduction in *x*<sup>10</sup> diameter is more pronounced at the higher shear velocity.

**Figure 8.** Influence of different shear velocities on the particle size density distribution.

#### **5. Conclusions**

In summary, it can be stated that crystal breakage always occurs in the case of combined normal and shear loading of lysozyme crystals. This is also true for normal stresses, for which no reduction in particle sizes can be detected without shear loading. In the case of combined stress, an increase in normal stress until a value of 64,288 Pa leads to significantly higher comminution. The shear velocity, on the other hand, has a much smaller effect on comminution and is a negligible influence compared to the normal stress. The high increase in fines can be explained by crystal abrasion in addition to crystal breakage. A

higher normal stress increases the friction in the sediment, which is why crystal corners and edges rub off. This contributes to the strong increase in fines.

Superposition of compressive and shear stress occurs in many technical apparatuses and must, therefore, be minimized in order to avoid undesirable comminution of protein crystals. Examples of an apparatus with superimposed stress include pumps, tubular centrifuges, decanter centrifuges, and cross-flow filters. With high crystal concentrations, such as those that occur during water evaporation crystallization, the risk of particle abrasion or particle collision is increased.

**Author Contributions:** Conceptualization, B.R. and H.N.; methodology, B.R.; validation, B.R.; formal analysis, B.R.; investigation, B.R.; resources, H.N.; data curation, B.R.; writing—original draft preparation, B.R.; writing—review and editing, B.R., M.G. and H.N.; visualization, B.R.; supervision, M.G. and H.N.; project administration, B.R. and H.N.; funding acquisition, H.N. All authors have read and agreed to the published version of the manuscript.

**Funding:** This research was funded by Deutsche Forschungsgemeinschaft grant number NI 414/26-2. We acknowledge support by the KIT-Publication Fund of the Karlsruhe Institute of Technology.

**Institutional Review Board Statement:** Not applicable.

**Informed Consent Statement:** Not applicable.

**Data Availability Statement:** The data presented in this study are available on request from the corresponding author.

**Conflicts of Interest:** The authors declare no conflict of interest.

#### **Abbreviations**

The following abbreviations are used in this manuscript:


#### **References**

