**5. Mitochondrial Function in Previously Frozen Specimens: Overview and Commentary on an Updated Approach to Respirometry**

Several attempts have been made to carry out respirometry on cryopreserved tissue samples. In these approaches, samples were collected and frozen in different solutions trying to preserve the structure and integrity of the mitochondrial membrane. However, the outcome of these attempts was variable and not standardized [143–147]. Mitochondrial complexes and supercomplexes can be isolated from frozen samples and separated in native gels. These isolated complexes and supercomplexes preserve their enzymatic activity and more interestingly, mitochondrial supercomplexes respire [45,46]. These works demonstrate that the electron transport chain is not destroyed by freeze-thawing and that previously frozen sample respirometry is feasible [46]. The challenge is to convert the assay performed in isolated mitochondrial supercomplexes to previously frozen samples, without the need to isolate the ETC components.

We have recently published a new approach to measure respirometry in frozen samples (RIFS) that reconstitutes maximal mitochondrial respiration in previously frozen samples [148,149]. Using RIFS, maximal oxygen consumption of the ETC is measured by providing physiological electron donors and acceptors with oxygen consumption as an integrated readout. It is important to consider that in previously frozen samples, the natural provider of electron equivalents, the tricarboxylic acid (TCA) cycle, is not directly feeding the ETC. For that reason, we need to bypass the components that are lost during the freeze-thaw step by providing electron donors directly to the ETC and correct for variable permeabilization of mitochondrial membranes, when necessary. For example, NADH is used rather than pyruvate + malate since NADH, a reducing equivalent that donates electrons to CI of the ETC, does not have a permeability barrier in frozen tissue. The protocol is straightforward in almost all cases except for tissues that contain a high proportion of fibers and collagen such as skeletal muscle or tissues that are rich in membranes such as the brain. In these tissues, an additional protocol step is required for either enzymatic digestion (skeletal muscle) or additional membrane permeabilization (brain) for the reagents and electron donors, particularly NADH, to be accessible to mitochondria.

The RIFS protocol was optimized and validated using the Seahorse XF96 Extracellular Flux Analyzer, but should be suitable for any respirometry platform. The RIFS assay includes injection of NADH to measure maximal CI respiratory capacity, succinate to measure maximal CII respiratory capacity, as well as TMPD to measure CIV activity. TMPD provides electrons to CIV by directly reducing cytochrome c. Inhibitors of CI, CIII, and CIV are used as controls in each assay to determine sensitivity to respiratory chain complexes. RIFS preserves 90–95% of the maximal respiratory capacity in frozen samples and can be applied to isolated mitochondria, cells, and tissue homogenates with high sensitivity. As an example in cells, Figure 2 compares fresh and frozen respiration in human lymphocytes from two donors. In fresh intact and permeabilized cells, these lymphocyte samples demonstrated low inter-individual variability and the expected response to injected compounds in all assays (Figure 2A,B). However, oxygen consumption rates were lower than expected in permeabilized cells possibly due to incomplete permeabilization. RIFS formulation allows for simultaneous measurements of respiration driven by electron entry through CI and CII while determining individual CIV activity. Additionally, CV, or ATP synthase, activity can be measured via ATP hydrolysis. In lymphocyte cell lysates, respiration driven through CI and CII, as well as CIV and CV activity can be measured with low inter-individual variability. Furthermore, CV is sensitive to dose-dependent inhibition with oligomycin (Figure 2C–E).

**Figure 2.** Extracellular flux analysis in fresh and frozen lymphocytes from two donors. Respirometry was measured in fresh (**A**,**B**) and frozen lymphocytes (**C**–**E**). (**A**) Oxygen consumption was measured in intact cells before and after injection of oligomycin (CV inhibitor), FCCP (chemical uncoupler), and antimycin A (CIII inhibitor) and rotenone (CI inhibitor). (**B**) Respiration in permeabilized cells with CII substrates. A final injection of TMPD was included to measure CIV. (**C**) RIFS protocol with succinate and rotenone to drive respiration through CII. (**D**) CV activity measured by media acidification as a result of ATP hydrolysis. (**E**) Oligomycin-sensitivity in the CV measurement. Abbreviations: AA: Antimycin A; ATP: Adenosine triphosphate; FCCP: Carbonyl cyanide p-trifluoromethoxyphenylhydrazone; oligo: Oligomycin; Rot: Rotenone; Succ: Succinate; TMPD: N, N, N', N'-tetramethyl-p-phenylenediamine.

To measure CV activity, the extracellular acidification (ECAR) is measured to monitor the acidification of the experimental medium in response to hydrolysis of injected, saturating levels of ATP, as previously described in Divakaruni et al. [150]. CV activity can be measured in previously isolated, frozen mitochondria as well as tissue and cell homogenates. In this application, the addition of the CV inhibitor oligomycin acts as an important control to determine that the increased ECAR after injection of ATP is a result of ATP hydrolysis and not non-specific acidification.

This ECAR approach will provide a measure of maximal hydrolytic capacity as a readout of CV enzymatic activity. However, it is important to note that ATP hydrolysis, both through ECAR measurements or with spectrophotometric assays, may underestimate CV deficiencies depending on the specific mutated subunit. For example, mutations in the mtDNA-encoded ATP synthase Fo subunit 6 that cause severe impairment in ATP synthesis may show little or no defects in ATP hydrolysis [151]. Therefore, while CV activity is a useful readout, it cannot replace measurement of ATP synthesis or respirometry in fresh samples or provide information on physiological levels of ATP hydrolysis.

The simplified sample preparation with RIFS, particularly with homogenates, and the 96-well format of the Seahorse XF96 Analyzer allows using significantly less biological material. RIFS does not require isolation of mitochondria or the use of detergents or plasma membrane permeabilizers, which simplifies the methodology and minimizes changes induced by partial or over permeabilization due to the use of detergents. Acin-Perez et al. provide validation in mouse, zebrafish, and human samples [148]. Another advantage of using total tissue lysates is that it accounts for tissue-specific mitochondrial function, where, for example, soleus showed significantly higher respiration than quadriceps, as previously described [152]. In homogenates and cell lysates, respirometry can also be normalized to

mitochondrial content to provide a measure of respiration per functional unit. In addition to traditional methods to measure mitochondrial mass, we have recently described a high throughput MitoTracker Deep Red protocol that requires minimal biological samples that can be run in parallel to RIFS and used to normalize respiration [148,149]. Since differences in mitochondrial content can account for differences in respiration rates in homogenates, the normalization step can differentiate between changes in mitochondrial content and function when homogenate samples or cell lysates are loaded per total protein. Normalization to mitochondrial content is not necessary with isolated mitochondria since respiration in these samples is already normalized to mitochondrial content when loaded per total protein.

The RIFS technique demonstrated that highly oxidative tissues such as heart, brown adipose tissue, and brain showed higher mitochondrial respiration rates per milligram of tissue whereas white adipose tissue, known for its low mitochondrial content, showed the lowest respiration rates. These results suggest that RIFS reveals physiological differences depending on the energy and metabolic demand of the tissue. When comparing RIFS to traditional methods used in frozen samples such as spectrophotometric enzymatic assays, the latter provides information about the maximal activity of the individual complexes or some of them combined, whereas RIFS allows for an integrative measure of the ETC from CI or CII to the natural electron acceptor that is oxygen. Additionally, RIFS often requires less sample to perform the assay compared to enzymatic assays and can show increased sensitivity to specific inhibition of the ETC complexes, as displayed in Figures 3 and 4 where results were obtained with half as much starting material compared with enzymatic assays. In addition, the assay is straightforward where the number of reagents needed is minimal in comparison to spectrophotometric assays. The simple, standardized protocol of RIFS could limit the inter-laboratory variability that has been reported with enzymatic spectrophotometric assays [153].

Snap freezing samples in liquid nitrogen is the best approach for RIFS. However, due to possible delays between sample collection, freezing, and storage for clinical samples, an initial concern was how quickly samples need to be frozen to preserve mitochondrial function. For that reason, different sampling procedures were tested to determine which approach best preserved maximal respiratory rates when samples could not be immediately flash-frozen in liquid nitrogen. Samples placed on ice for up to 3 h after collection as well as samples immediately stored at −20 ◦C preserved the integrity and function of mitochondria ETC complexes suggesting that this approach could be suitable for samples that require a delay in processing [148,149].

#### *5.1. Applications for Respirometry in Previously Frozen Samples*

There are numerous applications where RIFS can be used in preexisting, stored frozen samples. In adipose samples from pheochromocytoma patients, differences in maximal respiratory capacity that were observed in fresh clinical samples were maintained after long-term storage at −80 ◦C [154]. RIFS generated comparable data to those initially obtained in freshly isolated mitochondria including increased respiration through CI, CII, and CIV and increased mitochondrial content. Furthermore, RIFS made it possible to assay all the samples together compared to fresh samples that had to be run one at a time. This application could be particularly important in age-related studies since it enables a longitudinal design in which bioenergetic profiling can be monitored and compared over years across the same patients. RIFS could also be used as a functional test for patients with mitochondrial disease using less invasive samples. We tested the RIFS protocol in lymphoblastoid cell lines (LCLs), an immortalized B lymphocyte cell line, from patients with the mitochondrial disorder Leber optic atrophy (LHON). These cells were derived from patients with a mtDNA mutation in the NADH dehydrogenase subunit 4 gene. Both RIFS and spectrophotometric assays were able to detect impaired CI activity compared with control cells and revealed similar differences between the LHON lines (Figure 3A). Moreover, we were able to detect impaired NADH-driven respiration with half the amount

of sample or the same amount of sample while simultaneously getting information about CII-driven respiration and CIV activity, both of which were, for the most part, similar to the control cells with the exception of LHON LCL line 1 that showed some decrease in maximal CII respiration (Figure 3B). Comparing a citrate synthase assay to the MTDR assay to measure mitochondrial content, the results from both assays were comparable, but again RIFS required less starting material, fewer reagents, and was less time-sensitive than the citrate synthase assay (Figure 3C).

**Figure 3.** Comparison of RIFS and spectrophotometric approaches in Leber optic atrophy lymphoblastoid cell lines. (**A**) CI activity compared with the RIFS protocol. The spectrophotometric CI assay required 0.5 × 106 cells/well compared with 0.2 × 106 cells/well with RIFS. The data were also normalized to control cells as a more direct comparison of the two techniques. (**B**) CI, CII, and CIV results were obtained using RIFS. (**C**) Mitochondrial content measurements with citrate synthase activity and the MitoTracker Deep Red plate reader assay. Abbreviations: LCL: Lymphoblastoid cell line; LHON: Leber optic atrophy; RIFS: Respirometry in frozen samples.

**Figure 4.** Comparison of RIFS and spectrophotometric approaches to measure Complex I in human lymphocytes treated with rotenone prior to freezing. (**A**) CI activity compared with the RIFS protocol. The spectrophotometric CI assay required 1 × 106 cells/well compared with 0.4 × 106 cells/well with RIFS. The data were also normalized to control cells as a more direct comparison of the two techniques. (**B**) Compiled CI, CII, and CIV results from two donors were obtained using RIFS. Abbreviations: CI: Complex I; RIFS: Respirometry in frozen samples; ROT: Rotenone.

An additional potential application for RIFS is in monitoring environmental toxicants to which mitochondria are a key target [155]. For that reason, applications targeted to direct screening of compounds leading to mitochondrial toxicity can benefit from RIFS by testing the functional impact of environmental exposures within a large population [155,156]. We tested rotenone exposure in lymphocytes treated with different concentrations of rotenone before freezing and compared the RIFS protocol to a CI activity spectrophotometric assay. Both approaches were able to reveal impairment in CI, but the RIFS assay was more than twice as sensitive and required half as much starting material (Figure 4A). Additionally, given the high oxygen consumption rates with the RIFS protocol, it is likely that less material is required than was used in these proof-of-concept assays (0.4 × 106 cells/well). The requirement of less starting material could be particularly important when the amount of sample is limited. In addition to CI, RIFS also provide information on CII and CIV demonstrating that the rotenone-induced defect was limited to CI (Figure 4B).

#### *5.2. Limitations of RIFS*

While the RIFS protocol was not developed to replace traditional respirometry approaches, it is an alternative method when using fresh samples is not feasible due to sample size and availability, or immediate access to instrumentation is not possible, a problem particularly relevant for testing human samples. Table 5 provides a comparison of the advantages and disadvantages of measuring respiration in fresh and frozen tissue samples. Respirometry in fresh samples will undoubtedly provide the most complete information about mitochondrial function and the bioenergetics of a sample given the coupled nature of fresh samples, that is, respiration is coupled to ATP synthesis. However, RIFS can provide a convenient alternative when respiration in fresh samples is not feasible or as a complementary method to determine maximal respiratory capacity in an integrated, yet uncoupled, system.

**Table 5.** Advantages and disadvantages of respirometry measurements in fresh and frozen tissue samples. Abbreviations: ATP: Adenosine triphosphate; TCA: Tricarboxylic acid cycle.

Since freeze-thaw cycles will disrupt the integrity of the mitochondrial membranes and results in the loss of soluble matrix components, an important limitation is that coupled respiration or ATP synthesis cannot be measured with RIFS. These limitations should be considered when using this approach as it will not be sensitive to changes in substrate transport, substrate oxidation, ATP turnover, or TCA cycle flux. However, it has been demonstrated that disease-related phenotype often correlates with a decrease in maximal respiration [119,157–159] and that defects in coupled respiration are mostly related to some mechanism of uncommon drug toxicity [160]. Therefore, there are multiple applications for RIFS when testing fresh samples is logistically impracticable. While the data presented in this review are preliminary as proof of concepts for use of RIFS in circulating blood cells, the data suggest it could complement the current blood-based bioenergetics approaches. Furthermore, as the field continues to search for a biomarker of mitochondrial function, the RIFS approach provides an additional method for measuring maximal mitochondrial respiration and may be useful for non-invasive samples.

#### **6. Future Outlook: Respirometry in Non-Invasive Samples**

In the last decade, there have been substantial improvements in approaches to measure mitochondrial function in human samples. In particular, measurement of mitochondrial function in circulating platelets and leukocytes has elevated translational mitochondrial research and improved the prospects of clinical monitoring of mitochondrial bioenergetics. However, blood collection still needs to be completed in a medical setting and requires fresh processing. This processing can be intensive, especially when isolating sub-populations of PBMCs, which can present challenges in the ability to scale up for clinical trials. Additionally, blood-based bioenergetics do not always demonstrate changes in mitochondrial function observed in primary tissues [161], which suggests that subpopulation testing is an important consideration as subpopulation bioenergetic correlations may differ between tissues and diseases. More research is needed to understand the bioenergetic differences between subpopulations, population variability, and selection of the appropriate cell type so that this approach can be standardized. With the growing number of mitochondrial drug targets being developed and tested in clinical trials, there is a need for further minimally/non-invasive samples that can be scaled up for clinical testing and large-scale population studies. The RIFS protocol has the potential to be used for these applications if it can be optimized from non-invasive samples such as cells from urine or saliva, or buccal mucosa cells collected from cheek swabs, which can be easily obtained and immediately frozen without processing for later testing. Since RIFS requires minimal starting material, it could be modified for use in non-invasive samples that require

minimal processing. While it is still unknown whether RIFS will be a suitable approach for non-invasive samples, there is some evidence that buccal swabs may be an interesting foundation for the development of non-invasive samples.

Buccal mucosa cells represent a potential non-invasive sample that can be used for mitochondrial bioenergetics profiling. Buccal swab testing of electron transport chain activity is possible [162], but has been reported to show poor sensitivity and specificity for mitochondrial disease [41]. Frederiksen et al. showed that there is a significant correlation in mutation load of the common A3243G mtDNA point mutation in four different tissues: blood leucocytes, buccal cells, skeletal muscle cells, and urine epithelial cells [163]. They suggest that genetic testing of buccal cells and urine epithelial cells may provide an important supplement to the diagnosis of mitochondrial diseases, but also note that skeletal muscle was still the most informative tissue with the highest mutation load [163]. Buccal swabs from a 6-year old patient with suspected myoclonic epilepsy and ragged red fibers (MERRF) presented common mtDNA deletions and decreased ETC activities compared with controls [164]. In a study of 40 patients with suspected mitochondrial disease, when comparing enzymatic activities in buccal mucosa cells to muscle biopsies, there was a 77% correlation for an isolated CI deficiency, 86% for combined CI and CIV defects, and 100% in CIV. Overall, detection of mitochondrial impairment was 82% in buccal swabs compared to muscle biopsies [165]. Alterations in buccal mitochondrial DNA and function have also been reported in lung cancer patients and autism spectrum disorder patients (ASD) [166,167]. Additionally, there is evidence that buccal swab samples can be used to monitor mitochondrial enzymatic activity in response to treatment. In a study of ASD children with or without mitochondrial disease, CI, CIV, and citrate synthase activity in response to fatty acid and folate supplementation was measured in buccal extracts [168]. The results suggest that changes in CI and citrate synthase activities can be observed in buccal swabs in response to common interventions and that increased activity in these enzymes was more discernable in the mitochondrial disease subgroup.

While studies measuring mitochondrial function in buccal swabs are still limited, enough data demonstrate mitochondrial functional changes in buccal swabs to warrant further investigation. It is possible that an integrated readout, such as RIFS, would provide both better sensitivity and specificity for buccal samples. However, a potential limitation to buccal swab bioenergetics measurements could be bacterial contamination that decreases measurement specificity to the human samples [169,170]. Removal of bacteria from these samples will likely be an important step to improve sensitivity in the respirometry assay and ensure that outcomes are accurate. While this approach in buccal mucosa is still hypothetical and would require optimization and validation, it could provide a method to rapidly screen maximal mitochondrial respiratory capacity and provide a complementary approach to blood-based bioenergetic measurements and measurements in primary tissues.

#### **7. Summary**

Impaired mitochondrial function has been shown to contribute to the development and progression of various diseases. There is an unmet need to clinically monitor mitochondrial function. However, until now this has been limited due to the labor-intensive immediate processing of the samples. Blood-based bioenergetics have met some of this need but do not exclude the necessity of immediate processing. With the development of RIFS as a complementary technique to measure mitochondrial function in previously frozen samples, the possibility of clinical monitoring of mitochondrial function in samples collected at remote sites, or retrospectively in samples residing in tissue biobanks becomes more apparent. While these emerging approaches require further optimization, in some cases, and validation of population variability and sensitivity to detect differences between healthy and disease patients, they represent a potential step forward in translation bioenergetic monitoring in human samples.

**Author Contributions:** Conceptualization, L.S. and R.A.-P.; investigation, B.S., L.S. and R.A.-P.; writing—original draft, L.S. and R.A.-P.; writing—review and editing, B.S., C.B., L.S., R.A.-P. and O.S.S.; writing—revised draft, C.B., L.S. and R.A.-P. All authors have read and agreed to the published version of the manuscript.

**Funding:** O.S.S. is funded by American Diabetes Association 1-19-IBS-049.

**Institutional Review Board Statement:** Ethical review and approval were waived for this study due to the use of de-identified/commercially available samples.

**Informed Consent Statement:** Patient consent was waived due to the use of de-identified/commercially available samples.

**Data Availability Statement:** Data available upon request.

**Acknowledgments:** The authors would like to acknowledge Tom Grammatopoulos for providing materials. The authors would like to thank Corey Osto and Shai Rambod, as well as Ilan Benador, Daniel Dagan, Ajit Divakaruni, Marc Liesa-Roig, Marcus Oliveira, Martin Picard, and Amy Wang for scientific discussions and valuable advice.

**Conflicts of Interest:** The following competing interests did not influence the writing of this review: O.S.S. is a cofounder of Enspire Bio; O.S.S. and L.S. are cofounders of MyChon Health.

#### **Abbreviations**

Frequently Used Abbreviations: AA: Antimycin A; ADP: Adenosine diphosphate; Asc: Ascorbate; ATP: Adenosine triphosphate; CI: Complex I; CII: Complex II; CIII: Complex III; CIV: Complex IV; CV: Complex V; DH: Dehydrogenase; DNA: Deoxyribonucleic acid; ECAR: Extracellular acidification rate; ETC: Electron transport chain; FCCP: Carbonyl cyanide p-trifluoromethoxyphenylhydrazone; LCL: Lymphoblastoid cell line; LHON: Leber optic atrophy; NADH: Reduced nicotinamide adenine dinucleotide; nDNA: Nuclear DNA; mtDNA: Mitochondrial DNA; MTDR: MitoTracker Deep Red; OCR: Oxygen consumption rate; Oligo: Oligomycin; OXPHOS: Oxidative phosphorylation; PBMC: Peripheral blood mononuclear cell; RIFS: Respirometry in frozen samples; ROS: Reactive oxygen species; Rot: Rotenone; Succ: Succinate; TCA: Tricarboxylic acid cycle; TMPD: N, N, N',N' tetramethyl-p-phenylenediamine.

#### **References**


**Ion Exchange in Glass – The Changes of Glass Refraction**

Roman Rogoziński

154 Ion Exchange Technologies

*Chemistry,* Vol. 356, pp. 420-424.

*Total Environmental,* 293, 97-105.

*Rare Earths,* Vol. 28, pp. 111-115.

with thorin, *Talanta,* Vol. 49, pp. 929-935.

*Science and Technology*, Vol. 11, pp. 241-247.

Van der Walt T.N., Coetzee, P.P. (1996) Distribution coefficients and ion exchange behaviour of some elements with Purolite S-950 in hydrochloric acid, *Fresenius' Journal of Analytical*

Vigneau, O., Pinel, C., Lemaire, M. (2001) Ionic imprinted resins based on EDTA and DTPA derivatives for lanthanides(III) separation, Analytica Chimica Acta, Vol. 435, pp. 75-82. Vito, I.E.D., Masi, A.N., Olsina, R.A. (1999) Determination of trace rare earth elements by Xray fluorescence spectrometry after preconcentration on a new chelating resin loaded

Xu, X., Zhu, W., Wang, Z., Witamp, G.J. (2002) Distributions of rare earths and heavy metals in field-grown maize after application of rare earth-containing fertilizers, *The Science of*

Yan, Ch., Jia, J., Liao, Ch., Wu, S., Xu, G. (2006) Rare earth separation in China, *Tsinghua*

Yin, S., Wu, W., Zhang, B., Hang, F., Luo, Y., Li, S., Bian, X. (2010) Study on separation technology of Pr and Nd in D2EHPA-HCl-LA coordination extraction system*, Journal of*

Zhang, Z., Li, H., Guo, F., Meng, S., Li, D (2008) Synergistic extraction and recovery of cerium(IV) and fluorine from sulfuric solutions with Cyanex 923 and di-2-ethylhexyl

phosphoric acid, *Separation and Purification Technology,* Vol. 63, pp. 348-352

Additional information is available at the end of the chapter

http://dx.doi.org/10.5772/51427


**4.6 Main effect plot for bead width**

The main effect plot will help to determine the optimum value of the input parameter. Main effect plot is obtained for S/N ratio and mean for bead width by using Minitab 18 statistical software [8]. The main effect plot will represent significant the level of each input parameter as shown in **Figures 1** and **2**. Optimum value to obtained optimum bead width with their significant level is given in **Table 18**.

**Experiment no. Process parameter**

*DOI: http://dx.doi.org/10.5772/intechopen.93974*

**Table 12.**

**Signal-tonoise ratio**

Larger is better

Nominal is best

Smaller is better

**Table 13.** *S/N ratio.*

**175**

*Orthogonal array actual value.*

**The goal of the experiment Data**

Target the response and you want to base the signal-to-noise ratio on standard deviations only

Maximize the response Positive �10log 1

Minimize the response Non-negative with a

**characteristics**

Positive, zero, or negative

target value of zero

**Signal-to-noise ratio formulas**

*n*

i¼0̄yi

�10log.[P

�10log.1/n [ P i¼0yi2]

2/s2 ]

[ P *i*¼01*=y*i2]

**P B FT**

*Application of Taguchi Method in Optimization of Pulsed TIG Welding Process Parameter*


**4.2 Magnetic resonance imaging (MRI-brain)**

MRI-brain is more sensitive than CT-brain for detection of early hypoxic-ischemic brain injury. However, MRI-brain is not as convenient as CT-brain to obtain early in PCAS patients [46]. Diffusion-Weighted Imaging (DWI) sequences of MRI-brain are the most sensitive for cytotoxic injury from hypoxic-ischemic brain insults [47]. Restricted water molecules within ischemic brain tissue cause DWI restriction leading to hypersignal intensity appearance (**Figure 6**) [48]. DWI restriction threshold of 650 x 10−6 mm2 /s in more than 9 percent of brain volume determines poor outcomes [49]. Diffusion Tensor Imaging (DTI) plays a significant role in white matter tractography with a similar principle of intercellular water diffusion in DWI [50]. Fractional anisotropy, a DTI parameter, is a quantitative measurement for white matter abnormality [51]. Quantitative whole-brain white matter fractional anisotropy measured by DTI between days seven and 28 after cardiac arrest can predict long-term neurological outcomes [52, 53]. The fluid-attenuated inversion recovery (FLAIR) sequences of MRI-brain can also detect cytotoxic injury from hypoxic-ischemic brain insults [54]. The appearance of hypoxic-ischemic brain injury detected by FLAIR adds up the specificity to DWI in predicting unfavorable outcomes [55].

#### **Figure 6.**

*Imaging of a PCAS patient on five days after ROSC: CT-brain (A) showed no hypodensity lesion, DWI (B) showed hyperintensity in deep nuclei corresponded with hypointensity in Apparent Diffusion Coefficient sequences (C).*


*2.5. Intracellular ROS Assay*

HaCaT cells were seeded at a density of 1.5 × 105 in 8-well Lab Teck chamber containing DMEM supplemented with 10% FBS and were grown to 80% confluence. These cells were first treated with 1.5 μM Ectoine, followed by exposure to 3 J/cm2 UVA irradiation for the prescribed amount of time. Cells were washed with PBS and DCFH2-DA method was used to determine the intracellular ROS production using the Olympus Software solution software for each condition [21].


**4. Discussion**

The outbreak of the COVID-19 pandemic generated a global public health emergency [55–58]. Taylor et al. reported that psychological responses to previous epidemics and pandemics depended on individual vulnerabilities such as intolerance, insecurity, disease susceptibility, and anxiety [59]. Therefore, the authors of the present study conducted research on the assessment of the level of stress, the method of coping with stress, and the occurrence of symptoms of SS disorders in students of physiotherapy during the COVID-19 pandemic.

Multiple researchers emphasized the importance of mental stress in shaping the pathology of the SS [17–20]. Stressors with somatic consequences are largely dependent on the patient's personality type and the ability to cope with stress. The long-term effects of stressors breaks the body's adaptive mechanisms and lead to the accumulation of various disorders, including SS disorders [60]. According to the authors' research on a group of 188 physiotherapy students, the symptoms of SS disorders were found in 150 participants (80%). Research conducted by Emodi-Perlaman et al. showed that the COVID-19 pandemic had a significant adverse effect on the psychoemotional state of Israeli and Polish populations, leading to an increase in bruxism and TMD symptoms [61]. In the author's own research, symptoms of bruxism, i.e., tooth clenching, were found in 75 participants and teeth grinding in 34 (50% and 22.7%, respectively). In studies by Przyta ´nska et al., patients with high levels of stress reported parafunctional functions, i.e., awake bruxism, almost six times more often [62].

As it results from the conducted research, in the group of students with symptoms of SS disorders, higher scores were observed in the stress severity rating scale compared to the group without SS symptoms. The difference between the groups was statistically significant, which may indicate a significant role of stress in the development of SS disorders. De Medeiros et al. conducted research among Brazilian medical students and indicated that social isolation and stressful situations caused by the pandemic may increase the number of people with TMD symptoms [63].

It is generally accepted that mental stress causes increase in the muscular tension in different parts of the body [64,65]. The body posture adopted in stressful situations related to the so-called "fight or flight" reaction affects muscle contractions within the SS, leading to pain in the head and neck area. Therefore, muscle contractions in the body's upper quadrant may be a part of a stress-related defensive behavior [14–17]. In our study, the majority of physiotherapy students reported at least one symptom of SS and pain. The most common symptoms were headache (68%), neck and shoulder pain (58%), and teeth clenching (50%). Many clinical studies seem to confirm the relationship between the occurrence of SS disorders and severe stress, mostly in young people entering adulthood [66–69]. According to Quintiliani et al., as many as 89.4% of the respondents experienced an increase in perceived stress during the pandemic (66% reported moderate stress and 23.4% high stress) [70].

The authors of the study, by using the PSS 10 scale, found that the average stress intensity among physiotherapy students was in the range of six sten (19 raw points), indicating the average intensity of stress. Similar conclusions were drawn by other authors assessing the intensity of stress in emergency medical students, who also obtained an average value of six [70]. Importantly, in the present study, as many as 46.3% of the respondents obtained results in the range of 7–10 sten, which, according to the current interpretation of the research tool (PSS 10), qualifies them as high stress subjects. It is also worth noting that the maximum score on the PSS 10 scale in the group with symptoms of SSD was higher than in the control group.

The relationship between stress and neurological disorders, and tension headaches and migraines has been expensively described [71–74]. The results of our own research confirmed the existence of positive correlations between the intensity of stress and the symptoms of SS disorders, because with the increase in the intensity of the sense of stress, the intensity of headache and pain in the neck and shoulder girdle also increased. The

above fact should prompt clinicians to evaluate the structures of SS in patients who report the above-mentioned pain symptoms in their history.

According to the literature, SS disorders and higher levels of stress are more common in women than in men [75–78]. Wang et al. found three times higher levels of depression and health anxiety in women than in men during the COVID-19 pandemic [79]. Moreover, in the research by Liu et al., women were stronger predictors of PTSS symptoms after the pandemic [80]. Our results are in accordance with the abovementioned findings, as we discovered that women with a symptom of SS disorders were characterized by a significantly higher level of stress and a type D personality (*p* < 0.005) compared to men (*p* > 0.05).

The analyses carried out by the authors of the publication with the use of the Mini-COPE test were aimed at examining the activity that students undertake in a stressful situation. These actions are known as "coping methods". Research conducted by Babicka-Wirkus et al. on a group of 577 students from 17 Polish universities showed that, during the COVID-19 pandemic, students most often used coping strategies such as acceptance, planning, and seeking emotional support [81]. The analysis of the authors' results did not bring similar conclusions, because out of 14 strategies, four were significant predictors, i.e., positive thinking strategies, discharge, taking psychoactive substances, and blaming oneself. The results also showed that the higher the level of perceived stress, the more frequently the remaining three strategies were used. The level of perceived stress most strongly influenced the frequency of using the discharge strategy, which is a strategy focusing on revealing negative emotions and is associated with a feeling of mental discomfort. The results obtained by the authors allow for the conclusion that maladaptive stress coping strategies in students, especially during the pandemic, may have long-term consequences for their psychophysiological health and academic achievement. According to the research, the presence of SS symptoms did not contribute to changes in the way of coping as all coping strategies retained their statistical significance (*p* < 0.05).

As it was shown in the study, the intensity of type D personality manifestations was positively and moderately strongly associated with the general intensity of stress (rho = 0.690; *p* < 0.001). Students experiencing higher levels of stress exhibited more symptoms of type D personality, and those with more severe symptoms showed higher levels of stress. In a study by Cho et al. conducted among Korean students, it was observed that type D personality is related to the level of perceived stress. Additionally, it has been concluded that there are gender differences in type D personality, stress, and coping strategies [82].

By analyzing the available scientific literature and the results of our own research, it can be determined with high probability that the COVID-19 pandemic increases the intensity of stress and may influence the development of SS disorders.

More research is required to better understand the impact of COVID-19 stressors on stomatognathic disorders.

The major limitation of the present study is the fact that it was carried out in a very specific population. Therefore, the results cannot be extrapolated convincingly to the general population.


Preface

Chapter 7 **Desmoplastic Reaction in Biopsy Specimens of T1 Stage**

Chapter 8 **In vivo Optical Diagnosis of Polyp Histology: Can We Omit**

**Section 3 Complications of Colonoscopy: Minimizing the Risks of a**

Chapter 9 **The Major Complications of Colonoscopy: Sedation-Related,**

Muhammed Sherid, Salih Samo and Samian Sulaiman

**Section 4 Colonoscopy Is Not the End of the Way: Other Screening**

Robert J. Richards and Jerome Zhengrong Liang

**Hemorrhage Associated with Polypectomy and Colonic**

**Sub-Mucosal Invasion 161**

**VI** Contents

Marco Bustamante-Balén

**Screening Tool 195**

**Perforation 197** Paul Miskovitz

Chapter 10 **Complications of Colonoscopy 215**

**Technologies 241**

Chapter 12 **Virtual Colonoscopy 253**

Chapter 11 **Colon Capsule Endoscopy: Quo Vadis? 243** Samuel N. Adler and Cesare Hassan

**Colorectal Cancer Plays a Critical Role in Defining the Level of**

Shigeki Tomita, Kazuhito Ichikawa and Takahiro Fujimori

**Pathological Examination of Diminutive Polyps? 177**

Colorectal cancer (CRC) is the third more commonly diagnosed cancer in both sexes and it is also the third leading cause of cancer-related death among men and women in the United States.It represents a global public health problem and a heavy burdenfor any healthcare system. Fortu‐ nately, most CRCs originate from an adenoma, a precursor lesion that is easily accessible and re‐ movable. This is the rationale for CRC screening programs, which are aimed to diagnose CRC at an early stage or even better, to detect and resect the advanced adenoma before CRC has devel‐ oped. Recent evidence suggests that colonoscopy with polypectomy is efficient in reducing the in‐ cidence of and mortality from CRC. Moreover, the CRC resection rates have recently been shown to be reduced after the implementation of colonoscopy-based screening programs. As a conse‐ quence, CRC screening programs, either directly with colonoscopy or based on fecal occult blood test, are being implemented in most countries. The physicians who are devoted to endoscopy or who are involved in CRC screening programs face several challenges. They are requested to par‐ ticipate in the organization of population-based CRC screening programs; screening colonoscop‐ ies have to be performed accomplishing quality standards; difficult endoscopic resections are expected because performing colonoscopy in asymptomatic individuals makes possible an early diagnosis of advanced adenomas; and finally clinical decisions have to be made based on the in‐ formation provided by the number and histology of the detected adenomas.

In this background, this book providesaccurate and updated information on several topics to be faced when implementing a CRC screening program.It is divided in four sections, the first one dealing with the epidemiology of CRC and the rationale and challenges of implementing such a program both from a European and American point of view. The second section tackles the issue of managing the findings of a screening colonoscopy, commenting on the use of endoscopic sub‐ mucosal dissection to resect difficult polyps, the significance and management of the malignant polyp and the possibility of using image enhancing technologies to make a presumptive histologi‐ cal diagnosis of diminutive polyps reducing costs without losing efficacy. The third section sum‐ marizes the main risks of diagnostic and therapeutic colonoscopy and how to manage them.

The relevance of the topics covered in it and the updated evidence included by the authors turns this book into a very useful tool to introduce the reader in the amazing and evolving field of colonoscopy and CRC screening.

> **Marco Bustamante, MD, PhD** University Hospital La Fe Valencia, Spain

**Section 1**

**Overview of Colorectal Cancer. Implementation**

**of Colorectal Cancer Screening Programs**

**Overview of Colorectal Cancer. Implementation of Colorectal Cancer Screening Programs**

**Chapter 1**

**Colorectal Cancer**

http://dx.doi.org/10.5772/53524

**1. Introduction**

Kouklakis S. Georgios and Asimenia D. Bampali

Additional information is available at the end of the chapter

**1.1. Epidemiology – Clinical presentation-screening**

its development will be discussed here.

common in developed than developing countries.

Colorectal cancer (CRC) is a common and lethal disease. The risk of developing CRC is influenced by both environmental and genetic factors. Colorectal cancer is the third most common cancer worldwide. Clinical symptoms develop late in the course of the disease, and precursor lesions (adenomas) can be easily detected and removed. The disease is a candidate for early detection and prevention by screening. The epidemiology of CRC and risk factors for

**Epidemiology** — CRC incidence and mortality rates vary markedly around the world [1]. Globally, CRC is the third most commonly diagnosed cancer in males and the second in females, with over 1.2 million new cases and 608,700 deaths estimated to have occurred in 2008.Incidence and mortality rates are substantially higher in males than in females [2]. It is the fourth most common cause of cancer death after lung, stomach, and liver cancer. It is more

In the United States, both the incidence and mortality have been slowly but steadily decreasing. Annually approximately 143,460 new cases of large bowel cancer are diagnosed, of which 103,170 are colon and the remainder rectal cancers. Annually, approximately 51,690 Americans

**Incidence** — There is significant geographical variation in age-standardized and cumulative, 0-74 year incidence and mortality rates. Globally, the incidence of CRC varies over 10-fold. The highest incidence rates are in Australia and New Zealand, Europe and North America, and the lowest rates are found in Africa and South-Central Asia [5]. The highest incidence rate of CRC is estimated in the Czech Republic [39-42]. These geographic differences appear to be

and reproduction in any medium, provided the original work is properly cited.

© 2013 Georgios and Bampali; licensee InTech. This is an open access article distributed under the terms of the Creative Commons Attribution License (http://creativecommons.org/licenses/by/3.0), which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

© 2013 The Author(s). Licensee InTech. This chapter is distributed under the terms of the Creative Commons Attribution License http://creativecommons.org/licenses/by/3.0), which permits unrestricted use, distribution,

die of CRC, accounting for approximately 9 percent of all cancer deaths [6].


*Article* **Non-Stationary Contaminant Plumes in the Advective-Diffusive Regime**

**Iván Alhama 1, Gonzalo García-Ros 1,\* and Matteo Icardi 2**

**Abstract:** Porous media with low/moderate regional velocities can exhibit a complex dynamic of contamination plumes, in which advection and molecular diffusion are comparable. In this work, we present a two-dimensional scenario with a constant concentration source and impermeable upper and lower boundaries. In order to characterise the plume patterns, a detailed discriminated dimensionless technique is used to obtain the dimensionless groups that govern the problem: an aspect ratio of the domain including characteristic lengths, and two others relating time and the horizontal length of the spread of contamination. The monomials are related to each other to enable their dependences to be translated into a set of new universal abacuses. Extensive numerical simulations were carried out to check the monomials and to plot these type curves. The abacuses provide a tool to directly manage the contamination process, covering a wide spectrum of possible real cases. Among other applications of interest, they predict the maximum horizontal and transversal plume extensions and the time-spatial dependences of iso-concentration patterns according to the physical parameters of the problem.

**Keywords:** contamination plume; advection-diffusion; universal curves

### **1. Introduction**

There are many non-stationary scenarios in large extension water-saturated porous media in which the existence of both the regional velocity and molecular diffusion of a solute in the fluid are combined. The primary interest in eventual processes of pollutant spreading is to determine the spatial and temporal evolution of the contaminant plumes in order to plan actions of control. The environmental effects involve aquifer contamination, a problem that has been widely treated in the scientific literature since the middle of the last century. Specific aspects that range from the mathematical-physical theoretical description of the process [1] to real cases that generally cause socio-economic impacts [2] have been addressed.

The objective of this work, as in other works where similar methodologies have been applied [3–8], is to obtain dimensionless groups that govern the expansive dynamics of the plumes caused by the simultaneous effects of advection and diffusion. These groups will collect the preponderance of one phenomenon over another, and they are the monomials according to which the temporal evolution of the horizontal and transversal plume extensions can be described and represented graphically by means of a universal abacus.

From a practical point of view, these type curves or universal abacuses allows hydrogeologists and engineers to confront contamination problems in groundwater systems through direct analyses to gain rapid predictions using curve matching and interpolation techniques, making computer simulations or field testing unnecessary. Although the accuracy of the results depends, to a great extent, on the experimental parameters, which could be difficult to obtain when time or economic variables are present, the universal curves can

**Citation:** Alhama, I.; García-Ros, G.; Icardi, M. Non-Stationary Contaminant Plumes in the Advective-Diffusive Regime. *Mathematics* **2021**, *9*, 725. https:// doi.org/10.3390/math9070725

Academic Editors: Marco Pedroni and Lucas Jódar

Received: 4 February 2021 Accepted: 25 March 2021 Published: 27 March 2021

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

**Copyright:** © 2021 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/).

establish a range of spatial-temporal contamination extensions based on the reliability of the field data introduced as input, which must be set as a starting hypothesis. These type curves can also be used to obtain a preliminary idea of the spread of the contamination when planning in–situ trace tests, which require setting the location of control points [9,10] as well as delimiting the time-dependent perimeters when establishing different types of land use [11,12].

To limit the scope of the work, large 2D rectangular scenarios in which the extension of the contamination is far from reaching the boundaries of the domain have been devised. The regional flow has been implemented from left to right and the horizontal faces of the domain are impervious to both flow and solute, simulating a narrowly-confined layer of 1m in depth parallel to the horizontal surface. The methodology described in this work to obtain the proper monomials consisted of, firstly, applying the discriminated dimensionless technique to the governing equations to obtain the dimensionless groups (a standard objective of dimensional analysis [13–18]), secondly, defining the interdependences among monomials or unknowns using the Pi theorem [19] and, finally, providing extensive numerical simulations to verify the reliability of the groups and to depict the universal abacuses.

It is worth mentioning that in the first step of the protocol, the use of a standard non-discriminated technique to obtain dimensionless groups in the mathematical model that govern this problem is a topic of heated debate among researchers [20,21]. The classic techniques used to obtain the dimensionless groups, once the governing equation has been written in its dimensionless form, would not, in fact, lead to a proper characterisation of the problem. These techniques do not generally introduce the hidden magnitudes (which are the unknowns of interest) or the time factor into the process, so that such unknowns cannot be derived from the groups by applying the Pi theorem. Instead, the proposed dimensionless protocol (discriminated and normalised nondimensionalisation of the governing equation) leads directly to the least number of independent groups and to the function of dependence between the unknowns and the physical and geometric parameters of the problem. By means of mathematical approximations and manipulation, the protocol associates to each addend of the equation a numerical coefficient of a dimensional character which balances with the rest, deducing the groups as independent ratios between these coefficients.

It should also be noted that the dimensionless group that characterises this type of coupled problem (advection and molecular diffusion) is the so-called Peclet number [22]. A dimensional study of the equation has been carried out in some research with different transport and flow conditions in dispersive scenarios, depicting type curves that are dependent on the Peclet number [23,24]. Despite the fact that hydrodynamic dispersion effects can be neglected when regional velocities are small enough, in the scenario presented here, this number cannot be defined a priori since the extensions of the domain (length or width) are not relevant to the study of the dynamics of contamination. As mentioned above, in this work we consider instead the time needed for the contamination plume to reach the boundaries. The dimensionless groups deduced, all of which contain unknowns, allow us to establish the functional relationships that we are interested in.

As for the second step in the methodology, once the Pi theorem has been applied and the dependences of the unknowns have been deduced, the direct quotient between the regional velocity (vo) and the molecular diffusivity (D) emerges as a determining factor in the prevalence of the transport phenomenon over diffusion. It is a dimensional factor (vo/D) whose unit is the inverse of a length (m−1) with no apparent physical meaning. Indeed, the unknowns (horizontal and transversal extensions of the contamination plume or global isoline pattern) are also dependent on a time factor. The expression of this dependence allows us to separate the effect of the quotient (vo/D) and the product voτ (τ is time) to substantially simplify the universal representation.

Finally, the numerical simulations are performed with SEAWAT [25], a widely recognised and reliable software package used for theoretical [26,27] and practical-technical

purposes [28,29]. A large number of numerical simulations has allowed us to develop a wide set of universal abacuses that provide relevant information on the space-time dynamics of the contamination plume. The physical variables of the problem cover a range of solutions which are fully representative of all the cases that may occur in practice, including those asymptotic cases of negligible diffusion versus advection, and vice versa. The limit of negligible diffusivity allows us to represent new universal curves to characterise pure advective processes. The use of universal graphs is illustrated with examples.

The hypotheses assumed in the physical model are the following: (i) velocities small enough to neglect the effect of hydrodynamic dispersion, (ii) the Darcy flow (negligible inertial forces and laminar flow), (iii) isotropic and uniform hydraulic conductivity, (iv) viscosity independent of concentration, (v) two-dimensional geometry and the absence of gravitational effects, (vi) isothermal conditions, (vii) water-saturated porous media, (viii) the absence of sources and/or sinks for flow, (ix) fully miscible single-phase fluid with negligible compressibility for both the fluid and the porous matrix (constant porosity) and (x) non-reactive solute transport.

#### **2. Nomenclature**

#### **3. Mathematical Model**

Under the conditions mentioned at the end of Section 1, the mathematical model is composed of the mass conservation equations for the fluid and for the contaminant [30,31]. In mathematical terms, these read as follows:

$$-\left[\nabla \cdot (\rho \mathbf{q})\right] = \theta \frac{\partial(\rho)}{\partial \mathbf{\bar{r}}}\tag{1}$$

$$
\left[\nabla \cdot \mathbf{D}(\nabla \mathbf{c})\right] - \nabla \cdot (\mathbf{v}\mathbf{c}) = \frac{\partial \mathbf{c}}{\partial \tau} \tag{2}
$$

The relationship between the specific discharge and the actual fluid velocity is given by **q** = **v**θ, with θ being the porosity of the medium. Thanks to the coupling, according to Equation (2), the velocity of the fluid in the porous medium causes continuous redistribution of the concentration and, therefore, of the density of the solution, which affects its movement through Equation (1).

Equation (1) can be solved in terms of hydraulic potential or pressure by means of Darcy's law (i.e., a stationary momentum equation), expressed by Muskat [32] as follows:

$$\mathbf{q} = -\frac{\mathbf{k}}{\mu} \nabla \mathbf{p} \tag{3}$$

Since the global pressure distribution depends on the global density distribution, the qx and qy components of the specific discharge are affected by spatial variations in density. In terms of the hydraulic head, h = P ρg , and assuming that the fluid viscosity is not dependent on concentration, Darcy´s equation is written as **q** = −K*∂*h *∂*x , with K being hydraulic conductivity, a physical parameter of the porous medium, which is decisive for velocity, included in the mathematical model.

We consider a two-dimensional water-saturated rectangular domain with a 2:1 aspect ratio, according to the schematics of Figure 1.

**Figure 1.** Porous domain and boundary conditions.

The parameters of the soil that will later be shown to influence the values of the monomials have values of between 2.5 and 20 m/d (equivalent to fine sand or silty sand) for the hydraulic conductivity K, an effective porosity of 0.15 and a molecular diffusivity coefficient D ranging between 0.0003 and 0.0012 m2/d (representative of different salt species).

A zero-concentration regional flow, with a constant velocity vo and constant density ρo, enters from the boundary on the left due to a fixed hydraulic potential drop between the left and right-hand boundaries. The upper and lower horizontal boundaries are impermeable (zero normal velocity) and do not let the concentration go through by diffusion (zero

normal concentration gradient). In the position (x = xo, y = 0), which we will call focus, there is a cell of constant concentration co, representing the source of contaminant. The boundary conditions can be summarised as follows:

Hydraulic head:

$$\mathbf{h}\_{\text{(x=0,y,\tau)}} = \mathbf{h}\_{\text{l}\prime} \quad \mathbf{h}\_{\text{(x=L,y,\tau)}} = \mathbf{h}\_{\text{r}\prime} \quad \frac{\partial \mathbf{h}}{\partial \mathbf{y}\_{\text{(y=\pm\frac{M}{2},\tau)}}} = 0 \tag{4}$$

Concentration:

$$\mathbf{c}\_{(\mathbf{x}=\mathbf{x}\_0, \mathbf{y}=\mathbf{y}\_0, \mathbf{r}=0)} = \mathbf{c}\_0, \quad \mathbf{c}\_{(\mathbf{x}, \mathbf{y}, \mathbf{r}=0)} = 0 \quad \text{(at the rest of the domain)}\tag{5}$$

Under these conditions, we will study the non-stationary pattern of concentrations when H is large compared to the transversal dimension of the contamination plume. The isolines of the concentration pattern are deformed ellipsoids that widen in the horizontal and vertical directions of the domain due to the coupled processes of diffusion and advection. These ellipsoids surround the focus and finally keep its left end in a fixed position. As long as the contamination plume does not reach the right and upper (or lower) boundaries of the domain, the problem is not stationary (except for the small region to the left of the focus) so we expect the vertical and horizontal extent of the plume to depend on the parameters of the problem and on time. The interaction between the advection and diffusion phenomena is not easily predictable, even though each equation is linear in its unknown, as the position of the plume and isolines are complex non-linear functions of the velocities, pressure, and their gradients.

In the following section, we rewrite this mathematical model in its dimensionless form and, after some algebraic manipulations and the application of the Pi theorem, the dependences between the variables of interests and the physical parameters of the problem are obtained.

#### **4. Dimensionless Groups**

#### *4.1. Deduction of the Dimensionless Groups*

Before deriving the dimensionless groups for the coupled diffusion and advection problems or for purely advective problems, it is worth mentioning the easier problem of pure diffusion for which analytical solutions have been established [33]. Some comments should be made regarding the applied protocol.

The procedure followed to deduce the dimensionless groups that rule a given problem based on the Pi theorem consists of reducing the governing equations to their dimensionless normalised forms and obtaining such groups by comparing the coefficients that multiply the derivatives of each addend of the equation. These coefficients, which are physically or dimensionally homogeneous, are of the same order of magnitude in the normalised equation so their ratios are the dimensionless groups that are sought. In recent years, a formal procedure of nondimensionalising has been proposed and successfully applied to many complex problems in different areas to obtain new universal solutions. For example, it is worth citing the consolidation problem in soil mechanics, flow and (heat or solute) transport in porous media, and a variety of mechanical problems, all of which are coupled and nonlinear [7,8,34]. It should be noted that since the dimensionless numbers for most engineering or scientific problems are already found in the literature, the former procedure is generally avoided, and classical known numbers are directly used in the Pi theorem. This, however, generally leads to poorer predictive capabilities compared to the approach presented here.

The first step carried out to deduce the dimensionless governing equation is to choose a suitable list of reference quantities to define the dependent and independent dimensionless variables. These are chosen either from the input parameters of the problem or (as in this case) by introducing suitable unknowns whose order of magnitude will later be deduced as a consequence of the application of the Pi theorem. The only criterion with which these reference quantities should be chosen is that the range of values of the normalised variables is as bounded as possible to the interval [0, 1]. This criterion allows the derivatives of these variables to be averaged to a value of the order of the unit throughout the entire length of the scenario. In addition, the discrimination must be applied [35,36]. This means that vector variables such as position or velocity must be made dimensionless by means of different references according to their spatial direction. In this way, dimensionless groups commonly called "aspect ratios" or "form factors" do not necessarily emerge directly as independent groups.

As in other problems of similar complexity, the procedure requires good physical knowledge of the problem and some experience to find the unknowns introduced as references, their order of magnitude (by application of the Pi theorem) and their exact solution through the universal curves obtained by a large number of precise numerical simulations. This procedure has some advantages over the classical dimensional analysis. Firstly, the dimensionless groups emerging in the process, some containing unknowns and others without unknowns, are formally obtained and the relationship between them constitutes the direct application of the Pi theorem. Secondly, the groups have a unitary order of magnitude since they are obtained as balances between the addends of the governing equation. Groups of an order of magnitude higher (or lower) than unity, can be neglected in the governing equation, which simplifies the problem. Thirdly, this procedure incorporates the dimensionless physical parameters into the deduced monomials, thus reducing the global number of groups and making the characterisation more precise. The resulting groups are the proper parameters to represent universal solutions.

#### *4.2. Coupled Advection and Diffusion Case*

We assume isotropic molecular diffusivity, the porosity θ to be constant, and we neglect the transversal velocity (vy) and its spatial derivatives (since they are of an order of magnitude much lower than the regional velocity and its changes). Therefore, it follows that vx ≈ vo and *∂*vx *∂*x is negligible, and Equation (2) in rectangular 2D coordinates is reduced to:

$$\frac{\partial \mathbf{c}}{\partial \mathbf{r}} = -\mathbf{v}\_0 \frac{\partial \mathbf{c}}{\partial \mathbf{x}} + \mathbf{D} \frac{\partial^2 \mathbf{c}}{\partial \mathbf{x}^2} + \mathbf{D} \frac{\partial^2 \mathbf{c}}{\partial \mathbf{y}^2} \tag{6}$$

Dimensionless variables for x, y, c, vx, and τ are defined (discriminately) as

$$\mathbf{x'} = \frac{\mathbf{x}}{\mathbf{l\_x^\*}} \quad \mathbf{y'} = \frac{\mathbf{y}}{\mathbf{l\_y^\*}} \quad \mathbf{c'} = \frac{\mathbf{c}}{\mathbf{c\_0}} \quad \mathbf{r'} = \frac{\mathbf{r}}{\mathbf{r^\*}} \tag{7}$$

In these definitions, co and vo are known parameters, while l ∗ x, l ∗ y and τ∗ are unknown parameters related to each other. For a given time characteristic (τ∗), l ∗ x and l ∗ y define the extension of the solute plume from the constant concentration point. For example, they can be defined as the region in which the concentrations are above a certain percentage of co. With this choice, dimensionless variables may be considered normalised since they vary in the range [0, 1]. Substituting (7) in (6), the last equation becomes dimensionless:

$$\frac{\mathbf{c}\_{\rm o}}{\mathbf{r}^\*} \left( \frac{\partial \mathbf{c}'}{\partial \mathbf{r}'} \right) = -\frac{\mathbf{c}\_{\rm o} \mathbf{v}\_{\rm o}}{\mathbf{l}\_{\rm x}^\*} \left( \frac{\partial \mathbf{c}'}{\partial \mathbf{x}'} \right) + \frac{\mathbf{c}\_{\rm o} \mathbf{D}}{\mathbf{l}\_{\rm x}^{\*2}} \left( \frac{\partial^2 \mathbf{c}'}{\partial \mathbf{x}'^2} \right) + \frac{\mathbf{c}\_{\rm o} \mathbf{D}}{\mathbf{l}\_{\rm y}^{\*2}} \left( \frac{\partial^2 \mathbf{c}'}{\partial \mathbf{y}'^2} \right) \tag{8}$$

Assuming the derivatives between brackets to be of the order of one, four coefficients are found to describe the solution (or patterns) of the equation in the domain defined by l ∗ x, l ∗ y and time τ∗. These are:

$$\frac{\mathbf{c}\_o}{\mathbf{r}^\*}, \quad \frac{\mathbf{c}\_o \mathbf{v}\_o}{\mathbf{l}\_\mathbf{x}^\*}, \quad \frac{\mathbf{c}\_o \mathbf{D}}{\mathbf{l}\_\mathbf{x}^{\*2}} \quad \frac{\mathbf{c}\_o \mathbf{D}}{\mathbf{l}\_\mathbf{y}^{\*2}} \tag{9}$$

These coefficients have to be of the same order of magnitude since their terms in the equation balance each other. The independent ratios between these coefficients, chosen for suitability, are the dimensionless groups. We choose them as follows:

$$\pi\_1 = \frac{\frac{\text{c}\_0 \text{D}}{\text{l}^{\text{v}^2}}}{\frac{\text{c}\_0 \text{D}}{\text{l}^{\text{v}^2}}} = \frac{\text{l}\_\text{x}^2}{\text{l}\_\text{y}^{\*2}} \equiv \frac{\text{l}\_\text{x}^\*}{\text{l}\_\text{y}^\*} \quad \pi\_2 = \frac{\text{v}\_\text{o} \text{l}\_\text{x}^\*}{\text{D}} \quad \pi\_3 = \frac{\text{v}\_\text{o} \pi^\*}{\text{l}\_\text{x}^\*} \tag{10}$$

To characterise the solution in the better way, each of the unknowns l ∗ x, l ∗ y and τ∗ must appear only in one group, unless this is not possible (which is the case of l ∗ x). Therefore, group π3 can be substituted for the product of groups π2 and π3. This finally leads to the alternative solution:

$$
\pi\_1 = \frac{\mathbf{l\_y^\*}}{\mathbf{l\_x^\*}} \quad \pi\_2 = \frac{\mathbf{v\_o l\_x^\*}}{\mathbf{D}} \quad \pi\_3 = \frac{\mathbf{v\_o^2} \mathbf{r^\*}}{\mathbf{D}} \tag{11}
$$

Group π1 is an aspect ratio of the domain and provides information about the relation between l ∗ x and l ∗ y. The other groups allow us to relate time and the horizontal length of the transition or salt contaminated region measured from the focus. According to the Pi theorem, the solution π2 = Ψ(π3) leads to:

$$\mathbf{l}\_{\mathbf{x}}^{\*} = \frac{\mathbf{D}}{\mathbf{v}\_{\mathbf{o}}} \Psi \left( \frac{\mathbf{v}\_{\mathbf{o}}^{2} \mathbf{r}^{\*}}{\mathbf{D}} \right) \tag{12}$$

with Ψ an unknown function of its arguments. Writing the solution in the following way:

$$\mathbf{l}\_{\mathbf{x}}^{\*} = \left(\frac{\mathbf{D}}{\mathbf{v}\_{0}}\right) \Psi \left\{ \left(\frac{\mathbf{v}\_{0}}{\mathbf{D}}\right) \mathbf{v}\_{0} \mathbf{r}^{\*}\right\} \tag{13}$$

we obtain interesting and useful information. As expected, each iso-concentrated line of the pattern (defined by the dimensional value c ) depends on time (τ∗), regional velocity (vo) and molecular diffusivity (D), and for the same values of vo D and voτ∗, the concentration isoline is the same. Equation (13) shows detailed information about this kind of dependence. Firstly, keeping the ratio vo D constant, the patterns for each c are the same for all the times so that voτ∗ is also constant. This means that if we take two scenarios, the first with the pair of values (vo,D) and the second with the pair (aovo,aoD), the patterns of dimensionless isolines ( c ) at a given time τ∗ for the first scenario is the same as that of the second one at time τ∗ ao . In this way, scenarios of the same value of the ratio vo D , for a given time τ∗ contain the information from the patterns corresponding to all the regional velocity values. This allows us to depict a set of abacuses, one for each value of the ratio vo D , in which the extension of each concentration isoline, l ∗ x( c ) may be represented as a function of time τ∗, choosing c as the parameter of the abacus. To use this information, the time has to be suitably corrected according to the value of the regional velocity. In the following section related to the construction of the universal curves, this will be further clarified.

With regards to the characteristic transversal length, the group π1 = l ∗ y l ∗ x allows us to write l∗ y∼ l ∗ x, which is equivalent to the dependence

$$\mathbf{l}\_{\mathbf{y}}^{\*} = \left(\frac{\mathbf{D}}{\mathbf{v}\_{\mathbf{o}}}\right) \Psi \left\{ \left(\frac{\mathbf{v}\_{\mathbf{o}}}{\mathbf{D}}\right) \mathbf{v}\_{\mathbf{o}} \boldsymbol{\tau}^{\*} \right\} \tag{14}$$

with Ψ being a new function of its argument. The same discussion made above for l ∗ x applies here.

Finally, we can expect that the part of the pattern which is located to the left of the focus (xo, y = 0) will be a steady-state pattern after a relatively short time characteristic for which each isoline defined by c is located at a characteristic distance from (xo,0), l ∗ x(left). Indeed, deleting the addend *∂*c *∂*t from Equation (6) and assuming changes only in the x coordinate, the only emerging dimensionless group is π = vol ∗ x(left) D . Then, the solution is given as:

$$\mathbf{l}\_{\rm x(left)}^{\*} = \left(\frac{\mathbf{D}}{\mathbf{v}\_{\rm o}}\right) \tag{15}$$

The solutions provided by Equations (13)–(15), and the universal curves that make use of them, constitute an important management tool since they provide the level of global contamination of the domain and its gradation, at every instant.

#### *4.3. Pure Advection*

For this scenario, with the same assumptions as the coupled problem and using the same dimensionless variables defined in (7), Equation (2) is reduced to:

$$-\frac{\mathbf{c}\_{\rm o}\mathbf{v}\_{\rm o}}{\mathbf{l}\_{\rm x}^{\*}}\left(\frac{\partial \mathbf{c}^{\prime}}{\partial \mathbf{x}^{\prime}}\right) = \frac{\mathbf{c}\_{\rm o}\partial \mathbf{c}^{\prime}}{\mathbf{r}^{\*}\partial \mathbf{t}^{\prime}}\tag{16}$$

which leads to two coefficients, covo l ∗ x and co τ∗ and one single dimensionless group π1 = voτ∗ l ∗ x . This results in the following dependence:

$$\mathbf{l}\_{\mathbf{x}}^{\*} = \mathbf{v}\_{0} \mathbf{\tau}^{\*} \tag{17}$$

As in the former analysis, each isoline defined by its dimensionless concentration c has its particular solution. So, l ∗ x(c ) is the horizontal extension at time τ∗, going from the furthest point of the isoline to the constant concentration position imposed by the inner boundary condition. This length is proportional to the regional velocity but changes from one isoline to another

$$\frac{\mathbf{l}\_{\chi}^{\*}}{\mathbf{r}^{\*}} = \mathbf{v}\_{\circ} \,\,\Psi(\mathbf{c}') \tag{18}$$

The function Ψ(c ) is universal and may be depicted by a single numerical simulation. The region of concentration is a slender arrow in which the advancing fronts of the lower concentration isolines are ahead of the higher concentration fronts. The distance between the fronts of any pair of isolines c 1 and c 2, increases with time and depends on their concentrations according to the expression:

$$\mathbf{l}\_{\mathbf{x}}^{\*}\left(\mathbf{c}\_{1}^{\prime}\right) - \mathbf{l}\_{\mathbf{x}}^{\*}\left(\mathbf{c}\_{2}^{\prime}\right) \\ \quad = \mathbf{v}\_{\mathbf{o}} \pi^{\*} \left\{ \Psi\left(\mathbf{c}\_{1}^{\prime}\right) - \Psi\left(\mathbf{c}\_{2}^{\prime}\right) \right\} \tag{19}$$

This result can also be conveniently represented by a universal abacus.

#### *4.4. Perspectives*

It is interesting to discuss here what can be deduced from the previous results in finite domains and large time periods in which contamination reaches the boundaries of the scenario. Although this is not the main subject of this work, the treatment of finite scenarios adds new dimensionless groups and makes the solution more complex. For example, a finite scenario in the horizontal direction (but very large in the transversal direction) introduces the extension L in the dimensionalising process, giving rise to a new dimensionless group. In contrast, if the scenario is finite in the transversal direction (and very extensive in the horizontal one), the extension H of the domain must be considered, which also gives rise to the appearance of a new group. Finally, if the scenario is finite in both directions, the introduction of L and H would create two new dimensionless groups.

In this way, the introduction of any new condition like a sink or source (an injection or abstraction well), or reactive transport (degradation tracers), increases the number of monomials and, thus, the universal solutions are translated into a set of abacuses in a way that each one is set according to a specific value of a selected dimensionless group. A scenario similar to the one here presented but introducing an abstraction well located at a specific distance would be an interesting case for further research. From the practical point of view, this represents an often-used procedure in field testing to estimate the physical parameters of soils. The use of universal abacuses will contribute to obtaining these parameters by means of an inverse protocol.

#### **5. Verification by Numerical Simulations**

The numerical simulations have been carried out using the programme SEAWAT V.4 [37]. The following figures illustrate some of the simulated scenarios created to obtain the universal graphs. For example, Figure 2 shows the iso-concentration patterns in a large 2D scenario with a focus or constant concentration point located near the lefthand boundary, in which advection and diffusion effects are coupled. The data are: vo = 0.0006 m/d, D = 0.0006 m2/d and co = 1000 kg/m3. The concentrations of the isolines are 10, 200... kg/m3, which correspond to 10%, 20%... of co. The vertical extension of each pattern has been trimmed to reduce the size of the figure. As shown, the patterns progressively extend the contaminated region in vertical and horizontal directions. The deformation of the isolines (more pronounced for small concentration isolines) with respect to the circular shape that they would have with pure diffusion is due to the advective effect.

**Figure 2.** Iso-concentration patterns for different times in a large scenario with advection and diffusion. t = 500, 2000, 6000 and 10,000 days. The blue lines and grey network represent piezometric lines and cell extension (1 m2), respectively.

In contrast, the patterns for a scenario with only the advection effect is shown in Figure 3. The data are: vo = 0.05 m/d, D = 0 and co = 1000 kg/m3. In this scenario, the focus has been converted into a vertical segment to better appreciate the progress of the concentration fronts (of equal length to that of the segment) of each isoline. It can be seen that the velocities of each front depend on the concentration of the isoline, with values greater than the regional velocity (vo) for the small concentration isolines and lower for those with higher concentrations.

Figures 4 and 5 show the evolution of the isoline fronts in x direction (for vo = 0.05 m/d) when D = 0 and advection dominates. The first one shows the temporal evolution of the width of the region of variable concentration; that is, the area affected by significant contamination, defined as the distance between the isolines of concentration 100 and 900 (10% and 90% of co, respectively). The solution for an identical scenario with a 2co concentration has been superimposed. The separation between isolines of the same relative concentration (10% and 90%) has the same evolution. The different lines depicted in Figure 5 establish the distance to the focus of any iso-concentration line (defined as a percentage in respect to the contamination source) as a function of time. As previously mentioned, the velocity of each front is constant but increases as the concentration diminishes; that is, the variation in

the inclination of different lines indicates that low concentration isolines spread faster than high concentration values.

**Figure 3.** Iso-concentration patterns for times 50, 200, 500 and 1000 days in a large horizontal scenario due to advection (without diffusion).

**Figure 4.** Distance between the fronts of concentrations 0.1co and 0.9co as a function of time. co = 1000 and 2000 kg/m3, vo = 0.05 m/d.

**Figure 5.** Location of the isoline fronts as a function of time. vo = 0.05 m/d.

Continuing with the illustration of the coupled diffusion and advection effects, Figure 6a–c show the typical concentration profiles for different ratios of vo/D and co = 1000 kg/m3. These figures can be contemplated as a different and more complete configuration of Figures 2 and 3; the x-y axis distribution of concentration for specific times is now plotted on the vertical axis of concentration (z coordinate) in Figure 6a–c. Figure 6a is the case of no diffusion, with c = 1000 kg/m3 and vo = 0.0006 m/a, a velocity that corresponds to the front of the isoline of concentration 0.5co. The profiles gradually decrease their negative slope, increasing the distance from the small concentration isolines to those of higher concentration. Below the legend, the spreading of contamination in the x-y coordinates corresponding to the curve t = 6000 days is represented. Positions of concentrations 900 and 100 are represented as normalised with the 0.9 and 0.1 values in coordinate z (vertical axis in the chart). Figure 6b is a comparison between the profiles after 80 days for a pure diffusion process (D = 0.0006 m2/d), pure advection (vo = 0.0006 m/d) and coupled advection-diffusion (D = 0.0006 m2/d, vo = 0.0006 m/d). The profile corresponding to the combination of both effects presents two marked inflection points, depending on the vo/D ratio. This is a consequence of the coupling between both effects.

Finally, Figure 6c shows the profiles for a coupled problem with D = 0.0006 m2/d and vo = 0.0003 m/d for different times. For this ratio, vo/D = 2 as well as for ratios that are larger than unity, the profiles present an interesting result. If we consider the same concentration value in all the profiles (for example 0.3, which corresponds to c = 300 kg), the point of the domain with this concentration moves to the right more and more slowly until it stops at a certain distance from the focus at a certain time. This distance is set by the curve corresponding to the last time studied. The lower the chosen concentration, the greater the distance and time at which the concentration is fixed, and vice versa. As will be seen in the next section, similar phenomena occur in the transversal direction. As in Figure 6a, an x-y coordinate image of iso-concentration lines for t=14,000 days has been included below the legend so it can be compared with the concentration values of the same line in the z coordinate (or vertical axis in the chart).

**Figure 6.** Concentration profiles. (**a**) Only advection, c = 1000 kg/m3 and vo = 0.0006 m/d. (**b**) For τ = 80 days, pure diffusion (D = 0.0006 m2/d), pure advection (vo = 0.0006 m/d) and advection-diffusion (D = 0.0006 m2/d, vo = 0.0006 m/d). (**c**) Coupled problem, D = 0.0003 m2/d and vo = 0.0006 m/d.

All these results, which have been qualitatively described, will be represented by universal curves in the following section, after a large number of numerical simulations have been made. All this will verify the mathematical dependences derived through the non-dimensionalisation process followed in Section 4.

#### **6. Universal Solutions**

#### *6.1. Scenarios with Only Advective Flow*

For these scenarios, the only universal curve comes from expression (18), which we rewrite in the form of:

$$\frac{\mathbf{l\_{\chi}}^{\*}/\mathbf{\tau^{\*}}}{\mathbf{v\_{0}}} = \frac{\mathbf{v\_{c'}}}{\mathbf{v\_{0}}} = \Psi(\mathbf{c'}) \tag{20}$$

This relationship represents the ratio between the velocity of the front of the dimensionless concentration c , and the regional velocity vo. c is the ratio between the actual concentration of the front (c) and the constant concentration of the source. Figure 7 shows the universal dependence vc vo , or Ψ(c ), on c . The c = 0.5 front travels at the regional Darcy velocity. For a wide range of concentrations, c ∈ [0.15–0.85], the velocities of the fronts deviate very little from the value of vo (less than 10%). Only for very small or very large values of c , the velocities are somewhat higher (up to 1.3vo) or lower (0.7vo) than vo, respectively.

**Figure 7.** Universal dependence vc vo = Ψ(c ) (only advective effect).

#### *6.2. Scenarios with Coupled Diffusion and Advection*

The universal curves presented in this sub-section derive from the expressions (13)–(15). The first expresses the longitudinal extension measured from the source of contamination of each concentration isoline, the second, the maximum transversal extension, and the third, the horizontal extension of the polluted region to the left of the focus. Each isoline is characterised by its dimensionless concentration, taking that of the focus as a reference. The way in which these characteristic lengths depend on the physical parameters of the problem, vo and D, allows us to organise the universal curves in the form of an abacus, with a specific value for the relationship vo/D for each one. Every curve in Figure 8 represents the horizontal extension (vertical axis) of the isolines of dimensionless concentration c (l ∗ x(c )) against a time factor τc (horizontal axis), for vo D = 1. The simulations to determine

the function of the dependence (13) have been carried out for values vo = 0.0006 m/d and D = 0.0006 m2/d. In this way, the time factor (τc) is related to real time (τ∗) by

$$
\pi\_{\rm c} = \left(\frac{0.0006}{\mathbf{v}\_{\rm o}}\right) \mathbf{r}^\* \tag{21}
$$

Thus, for the same ratio ( vo D =1) and different values of vo and D, the extension of any isoline (l ∗ x(c )) associated with a real time τ∗ is obtained by entering the abscissa axis with the value τc given by expression (21). For greater detail, the time factor has been separated into two intervals with ranges [100–2000 days] and [10–100 days], Figure 8a,b, respectively. For longer times, the contaminated regions stabilise for successively increasing values of c . The lower the concentration is, the sooner the stabilisation occurs. Figure 8c,d show the extensions of some isolines for times of [2000–20,000 days] and [2000–8000 days], respectively.

The abacus corresponding to the monomials vo D = 2, 5 and 10 are shown in Figures 9–11, respectively, with details similar to the previous ones. Numerical solutions have been obtained by retaining the regional velocity (vo = 0.0006 m/d) and changing the diffusivity to successive values D = 3·10−4, 1.2·10−4 and 3·10−4 m2/d. In this way, the time factor is related to real time through the same expression as the former abacus (21).

**Figure 8.** *Cont.*

**Figure 8.** Universal curve l ∗ x(c ) as a function of time factor τc. Parameter of the abacus vo D = 1, vo = 0.0006 m/d. (**a**) τc ∈ [100–2000 days], (**b**) τc ∈ [10–100 days], (**c**) τc ∈ [2000–20,000 days], (**d**) τc ∈ [2000–8000 days].

**Figure 9.** *Cont.*

**Figure 9.** Universal curve l ∗ x(c ) as a function of time factor τc. Parameter of the abacus vo D = 2. vo = 0.0006 m/d. (**a**) τc ∈ [100–2000 days], (**b**) τc ∈ [10–100 days], (**c**) τc ∈ [2000–20,000 days], (**d**) τc ∈ [2000–8000 days].

**Figure 10.** Universal curve l ∗ x(c ) as a function of time factor τc. Parameter of the abacus vo D = 5. vo = 0.0006 m/d. (**a**) τc ∈ [100–2000 days], (**b**) τc ∈ [10–100 days], (**c**) τc ∈ [2000–8000 days], (**d**) τc ∈ [2000–20,000 days].

**Figure 11.** Universal curve l ∗ x(c ) as a function of time factor τc. Parameter of the abacus vo D = 10. vo = 0.0006 m/d. (**a**) τc ∈ [100–2000 days], (**b**) τc ∈ [10–100 days], (**c**) τc ∈ [2000–20,000 days].

Similarly, the abacus corresponding to values vo D = 0.5 and 0.1, made with vo = 0.0003 and 0.00006 m/d, are shown in Figures 12 and 13, respectively. Accordingly, the time factors and real time are related by:

$$
\pi\_{\rm c} = \left(\frac{0.0003}{\text{v}\_{\rm o}}\right) \mathbf{\tau}^\*, \quad \pi\_{\rm c} = \left(\frac{0.00006}{\text{v}\_{\rm o}}\right) \mathbf{\tau}^\* \tag{22}
$$

**Figure 12.** *Cont.*

**Figure 12.** Universal curve l ∗ x(c ) as a function of time factor τc. Parameter of the abacus vo D = 0.5. vo = 0.0003 m/d. (**a**) τc ∈ [100–2000 days], (**b**) τc ∈ [10–100 days], (**c**) τc ∈ [2000–10,000 days], (**d**) τc ∈ [2000–20,000 days].

**Figure 13.** *Cont.*

**Figure 13.** Universal curve l ∗ x(c ) as a function of time factor τc. Parameter of the abacus vo D = 0.1. vo = 0.00006 m/d. (**a**) τc ∈ [100–2000 days], (**b**) τc ∈ [10–100 days], (**c**) τc ∈ [2000–20,000 days], (**d**) τc ∈ [2000–20,000 days].

In relation to the largest transversal extension (l ∗ y(c =0.1)) given by dependence (17), the curves depicted in Figure 14 show this parameter (vertical axis) for isoline c = 0.1, which represents 10% of the concentration of the focus, as a function of time factor τc (horizontal axis) for vo D = 0.1, 0.5, 1, 2, 5 and 10. As in the former figures, the time factor scale and real time are related by expressions τc = 0.0006 vo τ∗ for the curves vo D = 1, 2, 5 and 10, τc = 0.0003 vo τ∗ for the curve vo D = 0.5 and τc = 0.00006 vo τ∗ for the curve vo D = 0.1. The lack of monotony in the slope of the curves is a consequence of the coupling between advection and diffusion, which occurs at different times according to the relative influence of one effect or the other. The patterns tend to stabilise when the advective and diffusive processes balance each other, while the preponderance of the diffusive process clearly establishes an unsteady pattern. There is a relationship between this maximum transversal extension (l ∗ y(c =0.1)) and the horizontal location (l ∗∗ x(c =0.1)) measured from the focus (f) at which such extension occurs. Figure 15 shows this dependence, in addition to that of the dimensionless concentration c = 0.1 and the same ratios of vo D .

**Figure 14.** Universal curve l∗ y(c =0.1) as a function of time factor τc.

**Figure 15.** Universal curve l∗∗ x(c =0.1) as a function of time factor τc.

To finish, the universal curve corresponding to the location of the stationary contamination fronts to the left of the focus, defined by the extension l ∗ x(left,c ), is shown in Figure 16 for the general case vo D = 1, in which advection and diffusion process are equally balanced. According to Equation (15), the curve is independent of the regional velocity as long as the ratio vo D remains constant. The period of time necessary to achieve this stationary distribution of concentration profiles, estimated from different simulations, is expressed as t = 1.8 vo and varies from 750 to 6000 days, depending on the regional velocity (vo = 0.0024 to 6000 m/d, respectively). Further research may aim to establish new curves as a function of the vo D ratio as well as to define the characteristic stabilisation time as a function of this ratio and the regional velocity.

**Figure 16.** Universal curve c as a function of longitudinal position for vo = 0.0006 and 0.0012 m/d.

#### **7. Conclusions and Final Comments**

The parametric dimensionless characterisation following a discriminated and normalised dimensionless protocol of the governing equations has allowed us to deduce precise information about the dynamics of contaminant plumes in extensive 2D scenarios with a constant concentration focus under the effects of advection and molecular diffusion. The proposed non-dimensionalised procedure as well as the application of the Pi theorem have resulted in accurate expressions of the unknown functions of interest (that is, the horizontal and transversal extensions of the plume) on the physical and geometrical parameters of the problem. In the most complex case, the lengths that define the pattern of the plume extension depend on three parameters: the dimensionless concentration of the isoline, the ratio between molecular diffusivity and regional velocity and a time corrected by the regional Darcy velocity. Based on these results and thanks to a sufficient number of numerical simulations, the functions can be universally represented by means of an easy-to-use abacus that facilitates the monitoring and management of the contamination plume in most real cases.

In the proposed protocol for the search of the dimensionless groups, normalisation makes it possible to approximate the changes or partial derivatives of the dependent variables to the unit when these are averaged over the entire domain of the problem, while discrimination prevents the emergence of monomials of the type of geometric shape factors that unnecessarily increase the number of dimensionless groups.

The coupling of the diffusive and advective flows is not intuitive since the cross values of the variables and their derivatives are combined in the governing equation. Small diffusivities versus advection (to the magnitude of unity, at the most) advance the dragging effect by redistributing the concentration in the posterior or central area of the plume before it affects the diffusive effect. However, if diffusion predominates over advection, it will produce a redistribution of the concentration by diffusion and, over this dampened field of concentrations, advection occurs.

The scenario addressed here is only a sample of the problems of contaminant flow and transport in porous media. First, the geometry of the scenario may be different or finite, including the effects associated with gravity flow. In addition, scenarios with a constant initial concentration (not maintained), with a non-constant concentration or with a contaminated fluid injection well could be tackled. Finally, the study of the dynamics of contaminants under the effects of advection and hydrodynamic dispersion (in general, with negligible molecular diffusivity) is another pending issue to characterise and for which to propose new universal solutions.

**Author Contributions:** Conceptualisation, I.A. and M.I.; methodology, I.A., G.G.-R. and M.I.; software, I.A.; validation, I.A. and M.I.; formal analysis, I.A. and G.G.-R.; investigation, I.A.; resources, I.A.; data curation, I.A.; writing—original draft preparation, I.A.; writing—review and editing, I.A., G.G.-R. and M.I.; visualisation, I.A. and G.G.-R.; supervision, I.A. All authors have read and agreed to the published version of the manuscript.

**Funding:** This research received no external funding.

**Data Availability Statement:** Data openly available in a public repository that issues datasets with DOIs.

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

#### **References**


*2.2. LiDAR Data Collection*

#### 2.2.1. UAV Laser Scanning Data

The UAV-LiDAR data of the 16 plots were captured on 4 and 5 June 2020, with a Genius UAV LiDAR developed by SureStar, Beijing. The UAV LiDAR area coverage and plot distribution are shown in Figure 2. Each plot at the two sites had a fixed size of 30 m × 30 m. The UAV-LiDAR system consisted of a DJI M200 platform, an R-Fans-16 laser scanner, and an advanced navigation spatial dual-coupled GNSS and IMU sensor. The R-Fans-16 scanner features 16 scan layers with a 30◦ vertical field of view (FOV). It operates at a wavelength of 905 nm and supports 320,000 points per second. The measurement range of R-Fans-16 is 200 m. The flying height for data acquisition was 40 m above ground level with a speed of 5 m/s, resulting in an average point density of approximately 430 pts/m2.