Md. Nurul Amin^*, Nazma Parveen, Md. Abdul Alim
Department of Mathematics, Bangladesh University of Engineering and Technology, Dhaka, Bangladesh.
DOI: 10.4236/ajcm.2025.154024   PDF    HTML   XML   150 Downloads   632 Views   Citations

Abstract

This study explores the role of electrode-induced temperature variations in plaque and blood flow dynamics within atherosclerotic radial arteries using computer simulations based on finite element analysis (FEA). The radial artery represents a valuable avenue for in vivo assessment of systemic atherosclerosis, providing diagnostic and predictive insights into coronary artery disease (CAD). In this model, a typical human arterial pressure is applied in the outlet, resulting in changes to blood velocity, pressure, and heat distribution. Initial inlet velocities and outlet pressures are applied using a time-dependent sinusoidal function that mimics pulsatile blood flow. The model integrates three key equations: the Navier-Stokes equations to describe blood velocity and pressure distribution, electric current equations to simulate heat generation, and the heat equation to evaluate temperature changes in the arterial wall. The simulation results were validated by comparing the velocity values with previously published data on radial and ulnar artery flow. For temperature validation, the simulated thermal distribution at the plaque region was found to be consistent with the reported ranges in the studies by Shiqing Zhao et al. The numerical results of this simulation revealed significant temperature localization near the plaque adjacent to the electrode. Due to the presence of the plaque and surrounding tissue, increase in temperature in electrode affects the blood flow, resulting in a decrease in blood velocity. The localized temperature behavior shows a rapid rise, followed by a peak, and then gradual stabilization. Blood velocity follows a pattern consistent with normal radial pulse propagation but decreases before the plaque and increases afterward. These findings contribute to a deeper understanding of the hemodynamic and thermal behavior in atherosclerotic arteries. The study also suggests potential future applications, such as employing AI-controlled micro devices to apply localized heat in semisolid plaque conditions for therapeutic purposes.

Keywords

Atherosclerosis, Electrode, Radial Artery, Temperature, Blood Flow

Share and Cite:

1. Introduction

Atherosclerosis, a progressive vascular disease involving the accumulation of lipids, fibrous elements, and inflammatory cells within arterial walls, remains a leading cause of cardiovascular morbidity and mortality worldwide. This challenge is particularly pronounced in low- and middle-income countries such as Bangladesh, where limited access to early diagnostics contributes to high rates of coronary artery disease (CAD) [1]. Among peripheral arteries, the radial artery has gained prominence in both diagnostic and therapeutic contexts due to its superficial location and relatively consistent anatomy. Studies have shown that structural changes in the radial artery, such as intima-media thickening and calcification, are indicative of systemic atherosclerosis and predictive of adverse cardiovascular outcomes [2] [3].

Over the past two decades, computational modeling has become a crucial tool in cardiovascular research for analyzing complex interactions between fluid flow, arterial wall mechanics, and pathological conditions such as plaque buildup. Various modeling strategies have been employed to investigate atherosclerotic hemodynamics, including computational fluid dynamics (CFD), fluid-structure interaction (FSI) models, peristaltic flow analysis, finite volume methods (FVM), and finite element methods (FEM). For instance, Taylor et al. [4] established early frameworks for FEM-based simulation of arterial flow, while Gijsen et al. [5] applied 3D imaging-guided simulations to assess plaque-induced changes in wall shear stress. Alagbe et al. [6] modeled how velocity profiles contribute to plaque growth in carotid arteries, and Hariri et al. [7] employed FSI methods to simulate pulsatile flow in stenosed vessels.

Thermal effects in vascular systems are another important area of study, particularly when associated with electrode-based interventions such as ablation and impedance mapping. The seminal work of Tungjitkusolmun et al. [8] introduced a FEM model for radiofrequency (RF) tumor ablation using electrode-induced Joule heating. Although initially developed for hepatic applications, this approach offers a valuable foundation for modeling localized heating in vascular systems, as heat alters blood viscosity, pressure gradients, and wall shear stress. Recent advancements have extended this thermal modeling to atherosclerotic conditions. Wang et al. [9] proposed a conformal heating strategy to target plaque regions, demonstrating the potential of thermal methods in vascular therapy, where Shiqing Zhao et al. [10] studied a coupled thermal-electrical-structural model for balloon-based thermoplasty treatment of atherosclerosis. Other studies have modeled thermal effects under Joule and magnetohydrodynamic (MHD) conditions using nanofluid-based blood representations [11] [12]. Asha and Sunitha [13] examined peristaltic flow with Joule heating in a non-uniform channel, while Misra and Sinha [14] analyzed thermal radiation and MHD effects in permeable capillaries. Such studies highlight the complex interplay between thermal gradients and vascular hemodynamics, yet limited work has been done in modeling these effects specifically within the radial artery under atherosclerotic conditions.

Therefore, the current study aims to fill this gap by investigating the effects of electrode-induced localized heating on blood flow dynamics in an atherosclerotic radial artery. Using a finite element method (FEM) framework, adapted from tumor ablation models, we simulate how thermal gradients influence flow velocity, pressure distribution, and the extent of thermally induced tissue damage. This approach integrates bioheat transfer equations, electric current modeling, and vascular geometry, providing a realistic and high-fidelity simulation environment. By advancing current understanding of thermal-fluid interactions in diseased arteries, this work lays the groundwork for AI-controlled microelectrode technologies capable of delivering precise, non-invasive thermal therapy for plaque regression. It also contributes to the broader field of patient-specific vascular modeling by combining mechanical and thermal insights within a single simulation framework.

2. Methods

2.1. Simulated Geometry

Figure 1(a) presents the cross-sectional view of an atherosclerotic radial artery, highlighting the anatomical features of the affected region. Figure 1(b) displays the simulated geometry, which was designed based on morphological data from Kim et al. [15] to replicate the artery’s structural and functional characteristics. Here arterial diameter is 2.7 mm. The computational domain is divided into four key regions: the blood, representing the flowing medium; the arterial wall, serving as the structural boundary; the electrode, used to generate localized heating for analysis; and the plaque, which serves as the primary site for examining the effects of atherosclerosis and thermal variations.

(a) (b)

Figure 1. (a) Radial artery atherosclerosis; (b) 3D geometrical model of simulation.

2.2. Governing Equations

Continuity equation:

$\frac{\partial u_x}{\partial x}+\frac{\partial u_y}{\partial y}+\frac{\partial u_z}{\partial z}=0$ (1)

The blood flow fluid dynamics were described using the Navier-Stokes equations and the equations for the conservation of mass. The fluid is assumed to be incompressible and homogeneous, with a Newtonian behavior. The blood was assumed to be Newtonian in this study, as its viscosity is almost constant in arteries; therefore, non-Newtonian effects were neglected. This simplification is commonly used for arterial flow because the shear rates in medium and large arteries are sufficiently high for blood to exhibit nearly constant viscosity. Although blood is inherently non-Newtonian. However, we acknowledge this as a model limitation and note that non-Newtonian effects may become more relevant in regions of low shear or at smaller vessel scales.

The Navier-Stokes equations can be expressed as follows:

X- Momentum equations-

$\rho \left(\frac{\partial u_x}{\partial t}+u_x\frac{\partial u_x}{\partial x}+u_y\frac{\partial u_x}{\partial y}+u_z\frac{\partial u_x}{\partial z}\right)=-\frac{\partial p}{\partial x}\text{}+\mu \left(\frac{{\partial }^2{u}_x}{\partial x^2}+\frac{{\partial }^2{u}_x}{\partial y^2}+\frac{{\partial }^2{u}_x}{\partial z^2}\right)+F_x$ (2)

Y- Momentum equations-

$\rho \left(\frac{\partial u_y}{\partial t}+u_x\frac{\partial u_y}{\partial x}+u_y\frac{\partial u_y}{\partial y}+u_z\frac{\partial u_y}{\partial z}\right)=-\frac{\partial p}{\partial y}+\mu \left(\frac{{\partial }^2{u}_y}{\partial x^2}+\frac{{\partial }^2{u}_y}{\partial y^2}+\frac{{\partial }^2{u}_y}{\partial z^2}\right)+F_y$ (3)

Z- Momentum equations-

$\rho \left(\frac{\partial u_z}{\partial t}+u_x\frac{\partial u_z}{\partial x}+u_y\frac{\partial u_z}{\partial y}+u_z\frac{\partial u_z}{\partial z}\right)=-\frac{\partial p}{\partial z}+\mu \left(\frac{{\partial }^2{u}_z}{\partial x^2}+\frac{{\partial }^2{u}_z}{\partial y^2}+\frac{{\partial }^2{u}_z}{\partial z^2}\right)+F_z$ (4)

where ρ is the blood density, $u_x$ , $u_y$ , $u_z$ are the blood velocity along x, y, and z axis, p is the blood pressure, µ is the dynamic viscosity of the blood, $F_x$ , $F_y$ , $F_z$ are body forces per unit volume. The fluid properties such as the density and the dynamic viscosity were assumed to be constant.

Their values were, respectively ${\rho }_{blood}=1060\text{kg}/{\text{m}}^{\text{3}}$ and ${\mu }_{blood}=0.0045\text{Pa}\cdot \text{s}$

Electric Currents (ec)

$\frac{\partial }{\partial x}\left(\sigma \frac{\partial V}{\partial x}\right)+\frac{\partial }{\partial y}\left(\sigma \frac{\partial V}{\partial y}\right)+\frac{\partial }{\partial z}\left(\sigma \frac{\partial V}{\partial z}\right)=\frac{\partial {\rho }_e}{\partial t}$ (5)

where, σ is the conductivity of the material, ${\rho }_e$ represents electric density.

V is the electric potential which initial value 22 V is used here. The 22 V electric potential was selected based on therapeutic voltage ranges commonly used as reported by Tungjitkusolmun et al. [8].

Heat transfer equations

$\rho c_p\frac{\partial T}{\partial t}+\rho c_p\left(u_x\frac{\partial T}{\partial x}+u_y\frac{\partial T}{\partial y}+u_z\frac{\partial T}{\partial z}\right)=k\left(\frac{{\partial }^{2}T}{\partial x^2}+\frac{{\partial }^{2}T}{\partial y^2}+\frac{{\partial }^{2}T}{\partial z^2}\right)+Q+Q_{bio}$ (6)

where ρ is the blood density, T is the temperature changes, k is the thermal constant, C_p denotes heat capacity of blood in constant pressure, Q total power dissipation density (Joule heating term) and $Q_{bio}$ is the metabolic heat source.

Heat capacity of blood in constant pressure C_p = 4180 J/(kg*K).

The electric potential equation (Equation (5)) was discretized using the Galerkin finite element method. The electric field-induced Joule heating term $Q=\sigma {\lfloor \nabla V\rfloor }^2$ was computed and coupled as a source term into the thermal energy equation (Equation (6)), allowing simulation of localized temperature elevation due to electrode-induced currents.

Transform into Finite Element Equations: Let us consider $N_{\eta }$ and $H_{\lambda }$ as weight function or linear shape function, equation with the finite element method. Two weight functions ($N_{\eta }$ for velocity components and temperature and $H_{\lambda }$ for pressure) are used due to the need for different interpolation orders to satisfy the Ladyzhenskaya-Babuska-Brezzi (LBB) condition. The weighted-integral technique of the governing equations is applied for deriving finite element equations as follows:

${\displaystyle {\int }_{\Omega }{N}_{\eta }\left(\frac{\partial u_x}{\partial t}+\frac{\partial u_y}{\partial t}+\frac{\partial u_z}{\partial t}\right)\text{d}\Omega }=0$ (7)

$\begin{array}{l}{\displaystyle {\int }_{\Omega }{N}_{\eta }\left(\frac{\partial u_x}{\partial t}+u_x\frac{\partial u_x}{\partial x}+u_y\frac{\partial u_x}{\partial y}+u_z\frac{\partial u_x}{\partial z}\right)\text{d}\Omega }\\ =-\frac{1}{\rho }{\displaystyle {\int }_{\Omega }{H}_{\lambda }\frac{\partial p}{\partial x}\text{d}\Omega }+\frac{\mu }{\rho }{\displaystyle {\int }_{\Omega }{N}_{\eta }\left(\frac{{\partial }^2{u}_x}{\partial x^2}+\frac{{\partial }^2{u}_x}{\partial y^2}+\frac{{\partial }^2{u}_x}{\partial z^2}\right)\text{d}\Omega }\end{array}$ (8)

$\begin{array}{l}{\displaystyle {\int }_{\Omega }{N}_{\eta }\left(\frac{\partial u_y}{\partial t}+u_x\frac{\partial u_y}{\partial x}+u_y\frac{\partial u_y}{\partial y}+u_z\frac{\partial u_y}{\partial z}\right)\text{d}\Omega }\\ =-\frac{1}{\rho }{\displaystyle {\int }_{\Omega }{H}_{\lambda }\frac{\partial p}{\partial y}\text{d}\Omega }+\frac{\mu }{\rho }{\displaystyle {\int }_{\Omega }{N}_{\eta }\left(\frac{{\partial }^2{u}_y}{\partial x^2}+\frac{{\partial }^2{u}_y}{\partial y^2}+\frac{{\partial }^2{u}_y}{\partial z^2}\right)\text{d}\Omega }\end{array}$ (9)

$\begin{array}{l}{\displaystyle {\int }_{\Omega }{N}_{\eta }\left(\frac{\partial u_z}{\partial t}+u_x\frac{\partial u_z}{\partial x}+u_y\frac{\partial u_z}{\partial y}+u_z\frac{\partial u_z}{\partial z}\right)\text{d}\Omega }\\ =-\frac{1}{\rho }{\displaystyle {\int }_{\Omega }{H}_{\lambda }\frac{\partial p}{\partial z}\text{d}\Omega }+\frac{\mu }{\rho }{\displaystyle {\int }_{\Omega }{N}_{\eta }\left(\frac{{\partial }^2{u}_z}{\partial x^2}+\frac{{\partial }^2{u}_z}{\partial y^2}+\frac{{\partial }^2{u}_z}{\partial z^2}\right)\text{d}\Omega }\end{array}$ (10)

$\begin{array}{l}{\displaystyle {\int }_{\Omega }{N}_{\eta }\left(\frac{\partial T}{\partial t}+u_x\frac{\partial T}{\partial x}+u_y\frac{\partial T}{\partial y}+u_z\frac{\partial T}{\partial z}\right)\text{d}\Omega }\\ =\frac{k}{\rho c_p}{\displaystyle {\int }_{\Omega }{N}_{\eta }\left(\frac{{\partial }^{2}T}{\partial x^2}+\frac{{\partial }^{2}T}{\partial y^2}+\frac{{\partial }^{2}T}{\partial z^2}\right)\text{d}\Omega }+\frac{1}{\rho c_p}{\displaystyle {\int }_{\Omega }{N}_{\eta }\left(Q+Q_{bio}\right)\text{d}\Omega }\end{array}$ (11)

For a scalar function $u$ the Laplacian ${\nabla }^{2}u$ can be rewritten using the Gauss theorem:

${\displaystyle {\int }_{\Omega }{N}_{\eta }{\nabla }^{2}u\text{d}\Omega }={\displaystyle {\int }_S{N}_{\eta }\nabla u\cdot n\text{d}S}-{\displaystyle {\int }_{\Omega }\nabla N_{\eta }\cdot \nabla u\text{d}\Omega }$ (12)

Applying Gauss’s divergence theorem to the second order derivative terms of the equations for generating the boundary integral terms associated with the surface tractions and heat flux. The equations become,

$\begin{array}{l}{\displaystyle {\int }_{\Omega }{N}_{\eta }\left(\frac{\partial u_x}{\partial t}+u_x\frac{\partial u_x}{\partial x}+u_y\frac{\partial u_x}{\partial y}+u_z\frac{\partial u_x}{\partial z}\right)\text{d}\Omega }+\frac{1}{\rho }{\displaystyle {\int }_{\Omega }{H}_{\lambda }\frac{\partial p}{\partial x}\text{d}\Omega }\\ -\frac{\mu }{\rho }{\displaystyle {\int }_{\Omega }\nabla N_{\eta }\cdot \nabla u_x\text{d}\Omega }=\frac{\mu }{\rho }{\displaystyle {\int }_S{N}_{\eta }{\tau }_{xx}\text{d}S}\end{array}$ (13)

$\begin{array}{l}{\displaystyle {\int }_{\Omega }{N}_{\eta }\left(\frac{\partial u_y}{\partial t}+u_x\frac{\partial u_y}{\partial x}+u_y\frac{\partial u_y}{\partial y}+u_z\frac{\partial u_y}{\partial z}\right)\text{d}\Omega }+\frac{1}{\rho }{\displaystyle {\int }_{\Omega }{H}_{\lambda }\frac{\partial p}{\partial y}\text{d}\Omega }\\ -\frac{\mu }{\rho }{\displaystyle {\int }_{\Omega }\nabla N_{\eta }\cdot \nabla u_y\text{d}\Omega }=\frac{\mu }{\rho }{\displaystyle {\int }_S{N}_{\eta }{\tau }_{yy}\text{d}S}\end{array}$ (14)

$\begin{array}{l}{\displaystyle {\int }_{\Omega }{N}_{\eta }\left(\frac{\partial u_z}{\partial t}+u_x\frac{\partial u_z}{\partial x}+u_y\frac{\partial u_z}{\partial y}+u_z\frac{\partial u_z}{\partial z}\right)\text{d}\Omega }+\frac{1}{\rho }{\displaystyle {\int }_{\Omega }{H}_{\lambda }\frac{\partial p}{\partial z}\text{d}\Omega }\\ -\frac{\mu }{\rho }{\displaystyle {\int }_{\Omega }\nabla N_{\eta }\cdot \nabla u_z\text{d}\Omega }=\frac{\mu }{\rho }{\displaystyle {\int }_S{N}_{\eta }{\tau }_{zz}\text{d}S}\end{array}$ (15)

$\begin{array}{l}{\displaystyle {\int }_{\Omega }{N}_{\eta }\left(\frac{\partial T}{\partial t}+u_x\frac{\partial T}{\partial x}+u_y\frac{\partial T}{\partial y}+u_z\frac{\partial T}{\partial z}\right)\text{d}\Omega }\\ -\frac{k}{\rho c_p}\left({\displaystyle {\int }_{\Omega }\nabla N_{\eta }\cdot \nabla T\text{d}\Omega }-\frac{1}{\rho c_p}{\displaystyle {\int }_{\Omega }{N}_{\eta }\left(Q+Q_{bio}\right)\text{d}\Omega }\right)=\frac{k}{\rho c_p}{\displaystyle {\int }_{S_w}{N}_{\eta }{q}_w\text{d}{S}_w}\end{array}$ (16)

Where the surface tractions (${\tau }_{xx},{\tau }_{yy},{\tau }_{zz}$ ), which represent the normal and tangential stress components along the outflow boundary S and the velocity components and fluid temperature or heat flux ($q_w$ ) that flows into or out from the domain along wall boundary $S_w$ .

All nodes are connected with both velocities and temperature. The corner nodes are also associated with pressure. This means that a lower order polynomial is chosen for pressure and which is satisfied through the continuity equation.

The velocity component and the temperature distributions and the linear interpolation for the pressure distribution according to their higher derivative orders in the differential equations as:

$u_x\left(x,y,z\right)=N_{\delta }{U}_{\delta }$

$u_y\left(x,y,z\right)=N_{\delta }{V}_{\delta }$

$u_z\left(x,y,z\right)=N_{\delta }{W}_{\delta }$

$p\left(x,y,z\right)=H_{\lambda }{P}_{\lambda }$

$T\left(x,y,z\right)=N_{\delta }{T}_{\delta }$

where $\delta =1,2,3,4$ is for velocity component and the temperature distributions; $\lambda =1,2,3$ is for pressure.

Now putting the element velocity component distributions, the temperature distributions, and pressure distribution equations the finite element equations can be written as follows:

$K_{\eta {\delta }^x}{U}_{\delta }+K_{\eta {\delta }^y}{V}_{\delta }=0$ (17)

$\left(K_{\eta \delta }{U}_{\delta }+K_{\eta \delta {\epsilon }^x}{U}_{\delta }{U}_{\epsilon }+K_{\eta \delta {\epsilon }^y}{V}_{\delta }{U}_{\epsilon }+K_{\eta \delta {\epsilon }^z}{W}_{\delta }{U}_{\epsilon }\right)+\frac{1}{\rho }{R}_{\delta {\epsilon }^x}{P}_{\epsilon }-\frac{\mu }{\rho }{K}_{\eta {\delta }^{xyz}}{U}_{\delta }=Q_{{\eta }^x}$ (18)

$\left(K_{\eta \delta }{V}_{\delta }+K_{\eta \delta {\epsilon }^x}{U}_{\delta }{V}_{\epsilon }+K_{\eta \delta {\epsilon }^y}{V}_{\delta }{V}_{\epsilon }+K_{\eta \delta {\epsilon }^z}{W}_{\delta }{V}_{\epsilon }\right)+\frac{1}{\rho }{R}_{\delta {\epsilon }^y}{P}_{\epsilon }-\frac{\mu }{\rho }{K}_{\eta {\delta }^{xyz}}{V}_{\delta }=Q_{{\eta }^y}$ (19)

$\left(K_{\eta \delta }{W}_{\delta }+K_{\eta \delta {\epsilon }^x}{U}_{\delta }{W}_{\epsilon }+K_{\eta \delta {\epsilon }^y}{V}_{\delta }{W}_{\epsilon }+K_{\eta \delta {\epsilon }^z}{W}_{\delta }{W}_{\epsilon }\right)+\frac{1}{\rho }{R}_{\delta {\epsilon }^z}{P}_{\epsilon }-\frac{\mu }{\rho }{K}_{\eta {\delta }^{xyz}}{W}_{\delta }=Q_{{\eta }^z}$ (20)

$\left(K_{\eta \delta }{T}_{\delta }+K_{\eta \delta {\epsilon }^x}{U}_{\delta }{T}_{\epsilon }+K_{\eta \delta {\epsilon }^y}{V}_{\delta }{T}_{\epsilon }+K_{\eta \delta {\epsilon }^z}{W}_{\delta }{T}_{\epsilon }\right)-\frac{k}{\rho c_p}{K}_{\eta {\delta }^{xyz}}{T}_{\delta }+\frac{k}{\rho c_p}{Q}_{\eta \delta }{T}_{\delta }=Q_{{\eta }^t}$ (21)

The coefficients in the element matrices are in the form of the integral over the element area and along the element edges S and $S_w$ as:

$K_{\eta \delta }={\displaystyle \int N_{\eta }{N}_{\delta }\text{d}\Omega }$ $K_{\eta {\delta }^x}={\displaystyle \int N_{\eta }\begin{array}{c}\frac{\partial N_{\delta }}{\partial x}\end{array}\text{d}\Omega }$

$K_{\eta {\delta }^y}={\displaystyle \int N_{\eta }\begin{array}{c}\frac{\partial N_{\delta }}{\partial y}\end{array}\text{d}\Omega }$ $K_{\eta {\delta }^z}={\displaystyle \int N_{\eta }\begin{array}{c}\frac{\partial N_{\delta }}{\partial z}\end{array}\text{d}\Omega }$

$K_{\eta {\delta }^{xyz}}={\displaystyle \int \nabla N_{\eta }\cdot \nabla u_z\text{d}\Omega }$ $K_{\eta \delta {\epsilon }^x}={\displaystyle \int N_{\eta }{N}_{\delta }\begin{array}{c}\frac{\partial N_{\epsilon }}{\partial x}\end{array}\text{d}\Omega }$

$K_{\eta \delta {\epsilon }^y}={\displaystyle \int N_{\eta }{N}_{\delta }\begin{array}{c}\frac{\partial N_{\epsilon }}{\partial y}\end{array}\text{d}\Omega }$ $K_{\eta \delta {\epsilon }^z}={\displaystyle \int N_{\eta }{N}_{\delta }\begin{array}{c}\frac{\partial N_{\epsilon }}{\partial z}\end{array}\text{d}\Omega }$

$R_{\lambda {\epsilon }^x}={\displaystyle \int N_{\lambda }\begin{array}{c}\frac{\partial N_{\epsilon }}{\partial x}\end{array}\text{d}\Omega }$ $R_{\lambda {\epsilon }^y}={\displaystyle \int N_{\lambda }\begin{array}{c}\frac{\partial N_{\epsilon }}{\partial y}\end{array}\text{d}\Omega }$

$R_{\lambda {\epsilon }^z}={\displaystyle \int N_{\lambda }\begin{array}{c}\frac{\partial N_{\epsilon }}{\partial z}\end{array}\text{d}\Omega }$ $Q_{{\eta }^{u_x}}=\frac{\mu }{\rho }{\displaystyle \int N_{\eta }{\tau }_{xx}\text{d}S}$

$Q_{{\eta }^{u_y}}=\frac{\mu }{\rho }{\displaystyle \int N_{\eta }{\tau }_{yy}\text{d}S}$ $Q_{{\eta }^{u_z}}=\frac{\mu }{\rho }{\displaystyle \int N_{\eta }{\tau }_{zz}\text{d}S}$

$Q_{{\eta }^w}=\frac{k}{\rho c_p}{\displaystyle \int N_{\eta }{q}_w\text{d}{S}_w}$ $Q_{\eta \delta }={\displaystyle {\int }_{\Omega }{N}_{\eta }{N}_{\eta }\left(Q+Q_{bio}\right)\text{d}\Omega }$

3. Initial and Boundary Conditions

The initial and boundary conditions used in this study were derived with reference to previously published simulations by Selmi et al. [16] and Choudhury et al. [17], with slight modifications to better match radial artery geometry and flow conditions. Specifically, inlet velocity and outlet pressure profiles were adjusted within physiological limits to reflect radial artery behavior more accurately. In this simulation, the inlet velocity of blood flow is set to 0.2 m/s, while the outlet pressure is maintained at 10,132 Pa multiplied by time-dependent function $f\left(t\right)$ to mimic pulsatile blood flow typical of the human cardiovascular system. Since the normal arterial pulse rate typically ranges between 70 and 100 beats per minute, an average pulse rate of 83 beats per minute is assumed. As a result, a time-dependent piecewise function $f\left(t\right)$ of two sub-domain to reflect pulsatile flow (1st sub-domain is for systole and 2nd sub-domain is for diastole) is applied over a period of 0.72 seconds [16].

$f\left(t\right)=\{\begin{array}{ll}0.72\times \sin \left(2.785\pi t\right),\hfill & 0\le t\le 0.35\hfill \\ 0.40\times \sin \left(2.22\pi \left(t-0.35\right)\right),\hfill & 0.35<t\le 0.72\hfill \end{array}$ (22)

A non-slip boundary condition is imposed on the arterial walls, assuming impermeability, thus preventing fluid penetration as it traverses the vessel. The initial temperature and the arterial wall temperature were set to 37˚C.

4. Numerical Methods

4.1. Materials

Table 1. Material properties of different parts of the model.

Property	Blood	Artery	Plaque	Electrode
Electrical conductivity [S/m]	0.667	0.333	0.36	5.76e7
Thermal conductivity [W/(m*K)]	0.543	0.512	0.63	386.47
Density [kg/m3]	1060	960	1220	8935
Heat capacity at constant pressure [J/(kg*K)]	4180	3600	3900	383.9
Dynamic viscosity [Pa*s]	0.0045	-	-	-

The simulation (Table 1) was performed in COMSOL Multiphysics 6.0 using a time-dependent solver. Time step: 0.01 s. Simulation duration: 0 - 10 s. Fully coupled solver was used with direct linear system solver (PARDISO). Material properties are chosen from Selmi et al. [16], Chowdhuri et al. [17], Changdar et al. [18], and Wang et al. [9].

4.2. Meshing

In this study, the finite element method was employed to solve the set of equations defining the problem. The geometric model has been discretized into tetrahedral elements, with a central region representing blood, as illustrated in the accompanying figure. The total count of elements in the mesh amounts to 204,189. The artery and thickness are clearly depicted in Figure 2 and Figure 3 shows grid independency test for pressure.

Figure 2. Geometry mesh.

Figure 3. Convergence of pressure for various elements size.

4.3. Validation

Direct validation data for the radial artery are very limited; therefore, an established aortic model is used only to verify the general flow behavior. The governing Navier–Stokes equations apply uniformly to both large and medium arteries, so matching aortic flow trends confirms that the numerical method is physically consistent.

Figure 4(a) and Figure 4(b) illustrate validation platforms between the results of present study and previous study to ensure the accuracy of the reported simulations when geometrical dimension was chosen from the former reports of Selmi et al. [16]. At 0.1 s maximum velocity found near the first branch of aorta and at 0.2s near the second branch. This shows the velocity profiles are nearly symmetrical. The validations revealed that the outcomes have a robust resemblance with one another. A close graphical agreement between the two profiles confirms that the current simulation accurately reproduces the established results, thereby validating the reliability of the present model.

(a)

(b)

Figure 4. (a) Velocity profile of aorta (Selmi et al. [16]); (b) Velocity profile of aorta (Present simulation).

5. Results and Discussion

(a)

(b)

Figure 5. (a) Point-wise temperature evolution at (−15, 0, −1); (b) Line profile of temperature along artery length from (−35, 0, −1) to (15, 0, −1).

Figure 5, point graph (a) illustrates that at t = 0, the initial temperature of 37˚C rises over time. At the 10-second mark, the highest temperature, reaching 46˚C, is observed at the point (−15, 0, −1). Figure 4, line graph (b) is depicted here to visualize temperature fluctuations within the plaque area. The graph spans from a point (−35, 0, −1) preceding the plaque to a point (15, 0, −1) following it. The temperature in the region prior to the plaque gradually increases, converging with the electrode area, and subsequently declines towards the area beyond the plaque.

Figure 6 shows the 3D graph on temperature changes around the area of plaque inserted electrode. The highest temperature found 46^oC within 10 s near the electrode. Isothermal area where plaque material is heated with electrode is shown in Figure 7 which is necessary to observe the changing condition of primary deposit in plaque area before calcification in the radial artery. The combined analysis of Figure 5 and Figure 6 demonstrates that temperature becomes highly localized within the plaque region due to the combined effects of the plaque composition and surrounding arterial tissues. This localized heating alters the thermal gradient near the stenosis site, influencing blood viscosity and leading to a corresponding reduction in flow velocity.

Figure 6. 3D graph on temperature changes.

Figure 7. Isothermal area around the electrode.

Figure 8 illustrates the effect of heat on blood flow. A noticeable decrease in velocity is observed with the application of heat, where the peak velocity reduces from 32.94 cm/s to 31.86 cm/s. This reduction can be attributed to the thermal influence on blood properties, as heat causes a slight decrease in blood density and an alteration in viscosity, leading to diminished flow velocity within the arterial lumen.

Figure 8. Velocity profile in atherosclerotic radial artery (without heat) and (with heat) at the point (−30, 1, 0.5).

Table 2 depicts the comparison between the present result of velocity in atherosclerotic radial artery with heat and previous published work Amin et al. [19]. This clearly shows that due to application of heat velocity decrease.

Table 2. Comparison of present data with previous published data in atherosclerotic radial artery.

Topics	Atherosclerotic Radial artery without heat [19]	Atherosclerotic Radial artery with heat (present data)
Vmax(cm/s) before Plaque at t = 0.15 s	21.59	21.24
Vmax(cm/s) after Plaque area at t = 0.0.15 s	22.65	21.59
Peak Velocity at the area of Plaque	34.94	31.86
Pmax(Pa) before plaque area at t = 0.15 s	7279.75	7280.10
Pmax(Pa) after Plaque area at t = 0.15 s	7187.30	7187.73
Peak Pressure at the area of Plaque	7467.33	7213.25

Figure 9 illustrates the effect of plaque formation on blood flow within the atherosclerotic radial artery. The blue line represents the atherosclerotic artery, while the red line corresponds to the normal artery. The peak velocity in the normal artery is observed to be 21.53 cm/s, whereas in the atherosclerotic artery it increases to 32.94 cm/s. This comparison clearly indicates that the presence of plaque causes a significant acceleration of blood flow within the stenosis region due to the reduced lumen area and resulting flow constriction.

Figure 9. Velocity profile in normal radial artery (without plaque) and atherosclerotic (with plaque) radial artery at the point (−30, 1, 0.5) with inlet velocity 20 cm/s.

Figure 10 illustrates blood velocity fluctuations along an artery over time (0 s to 0.6 s). Notably, peak velocity of 34.94 cm/s occurred at 0.2 s, while the lowest is 7.6 cm/s at 0.4 s. Evidently, the velocity displayed a rising trend followed by a decrease, succeeded by a minor increase within a narrower range, and concluding with a subsequent decline. These dynamics mirror the artery’s response to cardiac phases. The highest velocity is observed at the plaque area due to artery constriction caused by plaque, indicating the presence of atherosclerosis.

t = 0.1 s t = 0.2 s

t = 0.4 s t = 0.6 s

Figure 10. Velocity changes at different times.

The 2D graph (Figure 11) presents velocity variations at a point (−30, 1, 0.5) before and a point (1, 1, 0.5) after a plaque area. Notably, at 0.2 s, a prominent peak of 21.2410 cm/s in Figure 11(a) and 21.5926 cm/s in Figure 11(b) is observed. It shows the velocity slightly increase after plaque area. In the present study, the arterial wall was modeled as rigid because fluid–structure interaction (FSI) was not included. This simplification reduces computational cost but does not fully represent the natural compliance of the arterial tissue. If wall elasticity were considered, the lumen would deform slightly during pulsatile flow, which generally leads to lower peak velocities, reduced wall shear stress (WSS), and smoother pressure wave propagation. Therefore, the rigid-wall assumption may cause overestimation of peak velocity and WSS, especially near the stenotic region. Incorporating FSI in future studies would provide a more physiologically realistic representation of arterial behavior.

(a)

(b)

Figure 11. (a) 2D graph at a point (−30, 1, 0.5) before plaque area; (b) 2D graph at a point (1, 1, 0.5) after plaque area.

The provided Figure 12(a) and Figure 12(b) illustrate velocity changes along a specific arc length. Arc length is measured along the central axis of the curved artery from the inlet toward the outlet, crossing the plaque. While exhibiting a somewhat intricate pattern, it is evident that the velocity gradually diminishes from its peak at 34.9 cm/s, interspersed with multiple fluctuations. Along the x-axis, velocity fluctuations indicate the presence of plaque effects, while along the z-axis, velocity fluctuations reveal the highest velocity occurring at the midpoint of the artery due to the non-slip condition of blood within the artery.

(a)

(b)

Figure 12. (a) 2D graph of velocity changes along arc length shown in figure. (b) 2D graph of velocity changes along cross section of artery shown in figure.

The pressure of blood flow with time is depicted in Figure 13. By examining the function of time with inlet and outlet pressure, we observe a clear similarity in the graphs. The maximum pressure of 7280.10 Pa in graph (a) and 7187.73 Pa in graph (b) indicates the pressure difference. Graph (a) shows the pressure before the plaque, while graph (b) shows the pressure after the plaque, confirming a pressure drop across the stenosis. This reduction in pressure means the heart must pump more forcefully to maintain adequate blood flow, increasing cardiac workload in atherosclerotic conditions.

(a)

(b)

Figure 13. (a) Point graph before plaque area. (b) Point graph after plaque area.

Table 3. Comparison of present data with previously published work for velocity.

Topics	Present data	Amin et al. [19]	Chowdhuri et al. [17]
V (cm/s) at t = 0.25 s	21.24	21.59	19.00

Our simulation (Table 3) yielded a peak velocity of 21.24 cm/s at t = 0.25 s, which closely aligns with the result reported by Amin et al. [19], who obtained a velocity of 21.585 cm/s in a comparable radial artery model without an embedded electrode. This close match supports the accuracy of our velocity field and highlights the isolated influence of electrode-induced thermal effects on blood flow. Furthermore, our results are consistent with those of Choudhari et al. [17], who observed a velocity of 19.00 cm/s in the ulnar artery. Although anatomically distinct, the ulnar and radial arteries share similar diameters, wall structures, and flow regimes, which justifies the use of the ulnar artery model as a comparative reference. These anatomical and hemodynamic similarities strengthen the credibility of our velocity validation.

For thermal validation, the simulated temperature distribution—showing peak values between 42˚C and 46˚C near the plaque—was compared with the findings of Shiqing Zhao et al. [10]. The observed temperature range in our model falls within the known therapeutic window, confirming that the localized heating is realistic and within safe limits. This supports the reliability of the coupled electro-thermal simulation. Experimental validation of the full thermo-fluid model, including in vitro or Doppler-based assessments, is planned for future work. The nuanced differences observed between our results and those from previous models underscore the effects of both atherosclerosis and electrode-induced heating on arterial hemodynamics, offering new insights into disease-influenced blood flow behavior. This study uses a single set of material properties for the plaque. In reality, plaque composition can vary significantly—lipid-rich, fibrous, or calcified—and these differences may alter thermal conductivity and heat transfer behavior. Such variations could influence the temperature distribution and hemodynamic response observed in this study. Future work should investigate how different plaque compositions affect the thermal and flow characteristics.

6. Conclusion

This study applied finite element analysis to explore the temperature in influencing blood flow dynamics within an atherosclerotic radial artery. The simulation results showed that externally applied heat caused a rapid initial rise in temperature near the plaque region, followed by a gradual stabilization. This thermal response directly affected the hemodynamics: blood velocity increased slightly before the plaque, peaked within the plaque due to narrowing, and remained elevated downstream. In contrast, pressure dropped where velocity increased, highlighting the inverse relationship commonly seen in stenosed arteries. The findings demonstrate that temperature is not merely a passive factor but an active driver of hemodynamic changes in a diseased artery. By accurately capturing these thermo-fluid interactions, the model offers a valuable tool for investigating heat-based therapeutic strategies. Specifically, it supports the potential use of controlled temperature application—possibly via AI-assisted micro devices—for localized treatment of early-stage atherosclerosis. Future work will focus on validating the simulation results through experimental and clinical data, and on extending the model to include more complex physiological variables. Such efforts will be essential to ensuring the safety and effectiveness of temperature-guided interventions.

Conflicts of Interest

The authors declare no conflicts of interest regarding the publication of this paper.

References

[1]	Li, Z., Tang, Z., Wang, Y., Liu, Z., Wang, G., Zhang, L., et al. (2022) Assessment of Radial Artery Atherosclerosis in Acute Coronary Syndrome Patients: An in Vivo Study Using Optical Coherence Tomography. BMC Cardiovascular Disorders, 22, Article No. 120.[CrossRef] [PubMed]
[2]	Eklund, C., Omerovic, E., Haraldsson, I., Friberg, P. and Gan, L.M. (2014) Radial Artery Intima-Media Thickness Predicts Major Cardiovascular Events in Patients with Suspected Coronary Artery Disease. European Heart Journal—Cardiovascular Imaging, 15, 769-775.[CrossRef] [PubMed]
[3]	Gaudino, M., Tondi, P., Serricchio, M., Spatuzza, P., Santoliquido, A., Flora, R., et al. (2003) Atherosclerotic Involvement of the Radial Artery in Patients with Coronary Artery Disease and Its Relation with Midterm Radial Artery Graft Patency and Endothelial Function. The Journal of Thoracic and Cardiovascular Surgery, 126, 1968-1971.[CrossRef] [PubMed]
[4]	Taylor, C.A., Hughes, T.J.R. and Zarins, C.K. (1998) Finite Element Modeling of Blood Flow in Arteries. Computer Methods in Applied Mechanics and Engineering, 158, 155-196.[CrossRef]
[5]	Gijsen, F.J.H., Wentzel, J.J., Thury, A., Lamers, B., Schuurbiers, J.C.H., Serruys, P.W. and Slager, C.J. (2004) A New Imaging Technique to Study 3-D Plaque and Shear Stress Distribution in Human Coronary Artery Bifurcations in Vivo. Journal of Biomechanics, 37, 683-692.
[6]	Alagbe, E., Amoo, T., Oriaifo, O. and Ayeni, A. (2024) Computational Simulation of the Effects of Blood Flow Velocity on Atherosclerosis Progression in a Human Carotid Artery. F1000 Research, 13, Article 426.[CrossRef]
[7]	Hariri, S., Soltanian-Zadeh, H., Pourabdollah-Nejad, D.S. and Ghaznavi-Ghoushchi, M.B. (2018) Numerical Simulation of Pulsatile Non-Newtonian Blood Flow through Stenosed Arteries Using Fluid-Structure Interaction. Journal of Applied Fluid Mechanics, 11, 1211-1221.
[8]	Tungjitkusolmun, S., Staelin, S.T., Haemmerich, D., et al. (2002) Three-Dimensional Finite-Element Analyses for Radio-Frequency Hepatic Tumor Ablation. IEEE Transactions on Biomedical Engineering, 49, 3-9.[CrossRef] [PubMed]
[9]	Wang, H., Zhao, S., Zou, J. and Zhang, A. (2023) A New Conformal Penetrating Heating Strategy for Atherosclerotic Plaque. Bioengineering, 10, Article 162.[CrossRef] [PubMed]
[10]	Zhao, S.Q., Wang, H.Y., Zou, J.C. and Zhang, A.L. (2023) A Coupled Thermal-Electrical-Structural Model for Balloon-Based Thermoplasty Treatment of Atherosclerosis. International Journal of Hyperthermia, 40, Article 2122597.[CrossRef] [PubMed]
[11]	Amin, M.N. and Parveen, N. (2021) Free Convection Flow with Joule Heating, Heat Generation and Viscous Dissipation along a Vertical Wavy Surface. AIP Conference Proceedings, 2324, Article 050011.[CrossRef]
[12]	El-Shahed, M. (2003) Pulsatile Flow of Blood through a Stenosed Porous Medium under Periodic Body Acceleration. Applied Mathematics and Computation, 138, 479-488.[CrossRef]
[13]	Asha, S.K. and Sunitha, G. (2019) Effect of Joule Heating and MHD on Peristaltic Blood Flow of Eyring-Powell Nanofluid in a Non-Uniform Channel. Journal of Taibah University for Science, 13, 155-168.[CrossRef]
[14]	Misra, J.C. and Sinha, A. (2013) Effect of Thermal Radiation on MHD Flow of Blood and Heat Transfer in a Permeable Capillary in Stretching Motion. Heat and Mass Transfer, 49, 617-628.[CrossRef]
[15]	Kim, J.U., Lee, Y.J., Lee, J. and Kim, J.Y. (2015) Differences in the Properties of the Radial Artery Betweencun, Guan, Chi, and Nearby Segments Using Ultrasonographic Imaging: A Pilot Study on Arterial Depth, Diameter, and Blood Flow. Evidence-Based Complementary and Alternative Medicine, 2015, 1-7.[CrossRef] [PubMed]
[16]	Selmi, M., Belmabrouk, H. and Bajahzar, A. (2019) Numerical Study of the Blood Flow in a Deformable Human Aorta. Applied Sciences, 9, Article 1216.[CrossRef]
[17]	Choudhari, P. and Panse, M.S. (2016) Finite Element Modeling and Simulation of Arteries in the Human Arm to Study the Aortic Pulse Wave Propagation. Procedia Computer Science, 93, 721-727.[CrossRef]
[18]	Changdar, S. and De, S. (2015) Numerical Simulation of Nonlinear Pulsatile Newtonian Blood Flow through a Multiple Stenosed Artery. International Scholarly Research Notices, 2015, 1-10.[CrossRef] [PubMed]
[19]	Amin, M.N. and Parveen, N. (2024) Simulation Analysis of Blood Flow in Human Atherosclerotic Radial Artery. SSRN Electronic Journal.[CrossRef]