Radiofrequency ablation has been used for pain management since the mid-1970s [46] and has significantly evolved from a therapy that was mainly employed for mitigating neuropathic pain to one of the most promising and frequently applied techniques in clinical practices for alleviating axial spine and musculoskeletal pain [1,16]. Today, the application of RFA in treating pain has steadily expanded and accepted for the treatment of facet joint, sacroiliac joint dysfunction and osteoarthritis [8,9,10,11,12,13,15,16,17,18,19,20]. However, there are significant contradictions and inconsistencies in the reported clinical results on the efficacy of RF procedures. Most of the clinical studies available in the literature are retrospective studies or case studies limited to reporting high-quality randomized controlled trials [47,48,49,50,51,52,53,54,55]. Furthermore, inappropriate selection criteria and treatment parameters could result in poor treatment outcomes, whereas anatomical variations, which are still not well-established, could limit the accurate interpretation of the obtained results. Thus, additional large-scale clinical studies of RF application in pain management are needed with longer follow-up periods to demonstrate its efficacy along with quantification of associated long-term adverse effects.

There has been continuous development in RF delivery systems and protocols in the quest to increase the ablation volume, which will lead to enhanced efficacy of the RF procedure. Continuous RF is the conventional form of energy delivery, and it causes a decline in the transmission of pain signals by damaging both sensory and motor nerve fibers [1,19]. Three different types of electrodes are typically utilized for delivery of RF energy to the target neural tissue, viz., monopolar, bipolar and cooled RF electrodes (arranged in ascending order of ablation volume generation). Another mechanism of RF energy delivery to the target nerve is the pulsed RF procedure. As mentioned earlier, contrary to the continuous RF that causes tissue destruction by heat, the pulsed RF is virtually painless and does not lead to any neural coagulation and irreversible tissue damage. Instead, it leads to an alteration in the pain signal transduction of the nerve fibers [1,19]. Although the efficacy of pulsed RF has been well-documented, the exact mechanism of action in the mitigation of pain is not fully understood yet [14,18,54,56,57,58,59]. Pulsed RF has shown promise in treating neuropathic pain and certain other clinical conditions where continuous RF is potentially hazardous, such as radicular pain, headaches, chronic shoulder pain, knee and hip pain, axial low back pain, and peripheral nerve pain [18,56,57]. However, further clinical studies are required to quantify the exact mechanism responsible for pain mitigation and broaden the clinical utility of these interventional pain procedures.

In order to capture the nuances of the biological cells/tissues exposed to external forces, a multiscale modeling approach provides an efficient and cost-effective alternative for bridging the different scales during the computational analysis [60,61,62,63]. Thus, further refinements in the proposed pain management models can be attained by using multiscale modeling approaches, whereby the damage caused due to RF currents can be quantified at a cellular level. This can be accomplished by coupling predominant phenomena during RF procedures that occur at different scales. Such multiscale models would not only assist in a better understanding of the pre-existing bio-physical behavior during pain management therapy, but also help in predicting the mechanisms that remain elusive to date and in generating new hypotheses for quantifying small-scale effects. Furthermore, while dealing with the small-scale effects of biological cells exposed to RF currents, the development and usage of stochastic models is practically unavoidable [62,64,65]. Moreover, bones in general, and cortical bones in particular, are biological tissues considered as part of our computational domain in the present study. They may exhibit additional effects, based on coupled phenomena, that would be useful to incorporate in further developments of the presented model. Among these, piezoelectricity plays a special role (e.g., [66]). For example, for cortical bones, piezoelectric properties are often responsible for the coupling between macroscopic and micro/nanoscopic scales, which may provide additional insight into certain dysfunctions and diseases [66]. Such properties also provide a foundation for wider usage of these biomaterials in tissue engineering [67]. In describing piezoelectricity, we couple electrical and mechanical fields. The well-posedness of such models of coupled piezoelectricity, along with rigorous energy bounds, were derived by one of us in a series of earlier papers, e.g., [68]. This was done for the first time in a general, dynamic setting through the application of the Faedo–Galerkin procedure and generalized solution technique. Coupled electro-mechanical models have been developed and used in a wide range of applications [69,70,71,72,73,74,75,76,77].

Thermal field treatment requires special attention for problems like those considered in this study. Notably, the thermal effects in this study were quantified considering the classical Fourier’s law of heat conduction with the Pennes bioheat equation as presented in Equation (5), which assumes the infinite speed of heat propagation. However, in biological tissues that have non-homogenous inner structures, accounting for a finite speed of thermal disturbances becomes important, suggesting the existence of non-Fourier conduction with a delay ranging from 10 to 20 s [38,78]. Several studies have been reported providing a way to incorporate such non-Fourier heat transfer behaviors in their computational models, e.g., [23,38,79,80,81,82]. Moreover, the attainment of elevated temperatures during RF procedures can also result in thermo-elastic deformation, including thermal expansion and tissue shrinkage, which is interlinked with many complex small-scale effects, such as protein denaturation, collagen contraction and dehydration. Again, the exact mechanism of such associated effects at an elevated temperature within the biological tissue during thermal therapies are not completely elucidated yet, but significant recent developments have been devoted to this area of research utilizing both experimental and computational studies [23,78]. From a computational perspective, the coupling between thermal and mechanical fields, e.g., for elastic tissues such as muscles, etc., can be done by the development of coupled models of thermoelasticity, as well as efficient numerical methods for their solution, e.g., [83,84,85,86,87,88,89,90,91,92,93,94,95]. Moreover, the development of such models also includes complex nonlinear cases where numerous advances have been made in the improvement of numerical methodologies, e.g., [96,97,98,99]. The coupling of the thermoelasticity and piezoelectric model, as is the case for bone tissues, can be done by the development of piezothermoeleastic models [100,101,102]. This could lead to the development of fully coupled thermo-electro-mechanical models of thermal therapies [23]. Also, the development of such models is relevant to other areas of application, as well as in the development of new methods [103,104]. Furthermore, the proposed model of RF application for mitigating pain presented in this study assumed the quasi-static approximation of Maxwell’s equations (see Equation (1)), whereby the extent of variation of electric and magnetic fields is negligible. Consequently, the electromagnetic field is modeled only by considering the electric field component (neglecting the magnetic field effects) because the wavelength of the electromagnetic field at RF frequency of around 500 kHz is several orders of magnitude larger than the size of the active electrode. However, models exist in the literature that consistently include the magnetic field effects, e.g., [103]. The electrode-tissue and plastic-tissue interfaces (as presented in Figure 1) can be more rigorously modeled by incorporating electrode-tissue contact force [60,105,106,107,108,109], as well as other non-trivial effects, e.g., thermal degradation, spiking, etc., in polymeric materials [110,111,112].

Several studies have been reported in the literature for modeling blood perfusion within biological tissues, at both micro- and macro-vascular levels [78,113,114,115,116,117,118,119]. Notably, micro-vascular perfusion refers to the perfusion at a capillary (or small-scale) level while macro-vascular perfusion is associated with the heat-sink effect caused by large blood vessels [120]. The blood flow at a micro-vascular level within the biological tissue is typically modeled by utilizing the porous media theory, whereby the tissue is assumed to be comprised of two phases: the solid phase comprising of cells and the extracellular space, and the fluid phase comprising of capillary size blood vessels [94,95,118,119,121,122,123]. The blood flow within the large blood vessels ( 2 mm in diameter) is modeled by additional coupling of the fluid flow model with the proposed thermo-electric model presented in this study, whereby the geometry of the blood vessel within the computational domain can be incorporated either by including a cylinder or a vascular tree [113,114,115,116,117,124,125]. It is expected that further refinement of the model can be done by deriving the computational domain from actual patient-specific data, which will provide more rigorous analysis and would help medical practitioners to obtain more accurate and precise predictions of the treatment outcomes during the RF application in pain management. Thus, the coupled multiscale framework could assist us in quantifying the unknown biologically-relevant phenomena occurring at cellular and sub-cellular scales and lead to a better understanding of the associated intricacies of RF application in pain management.

The application of machine learning has been growing rapidly in the biological, biomedical and behavioral sciences. Importantly, both machine learning and multiscale modeling complement each other in creating more robust predictive models in the current field of research [126,127]. Recently, several studies have been reported in the literature that have explored the application of AI and machine learning algorithms in the field of thermal therapies [40,128,129,130,131,132,133,134,135,136]. The integration of machine learning with the coupled models could play a vital role in decision-making processes and the treatment planning stage of such procedures, e.g., by providing a priori information about electrode placement for enhancing treatment efficacy or by the real-time monitoring of the damage to the target tissue and other critical structures. Furthermore, a general framework of incorporating human factors into mathematical models of complex systems with control has been provided in [137,138]. This can be useful in the context of AI and the machine learning algorithms mentioned earlier in Section 2. Moreover, there has been a surge in the development of neural tissue models for capturing the transduction, transmission, perception and modulation of pain at molecular, cellular and neuron networks levels [78,139,140,141,142,143]. The aforementioned coupled multiscale thermo-electro-mechanical model can be readily integrated with the Hodgkin–Huxley neural model for predicting the treatment outcomes in terms of decline in the actual pain signals that can be coupled with the damage model presented in Equation (7). Such coupling of the neuronal models with the proposed model would assist in our better understanding of the molecular changes affecting the neuronal behavior, such as in quantifying the exact damage to the axons during the application of RF procedures for treating pain. Future studies can be conducted by incorporating the actual physiological neuronal geometries and modeling of biophysical phenomena at sub-cellular scale, viz., accounting for changes in the concentrations of potassium, sodium, calcium and magnesium at the membrane layer [78,144,145,146,147,148]. Such multiscale, multiphysics and fully coupled models will provide a better understanding of the molecular changes affecting the neuronal behavior, along with quantification of the mitigation of actual pain signals during RF procedures.