Figure 2 shows an overview of the RFCA for arrhythmias, including blood in the cardiac chamber and the porous cardiac tissue with the inserted (pushed) radiofrequency ablation catheter. The ablative catheter consists of the catheter body and the active electrode made from polyethylene (PE) and platinum-iridium alloy (Pt-Ir). Figure 3 shows the 3-D computational schematic with domains and boundary conditions. Table 1 lists the model thermal, electrical and mechanical properties, along with other parameters. The mechanistic model with a porous media approach for RFCA, combined with an analysis of electrical, thermal, poromechanical deformation, and fluid transport, is simulated to obtain the results in this study.

Fig. 2

A schematic of the 3-D Radiofrequency Cardiac Ablation (RFCA) procedure and the computational domain

Fig. 3

A detailed schematic of the RFCA computational model geometry and boundary conditions used for electrical, thermal, fluid flow, and deformation analysis

Table 1 Input parameters including electrical, thermophysical, mechanical, and other properties used in numerical computations of the RFCA model

List of assumptions made to simplify the model include: 1) there is no chemical reaction in the tissue, 2) the steam pop is described by tissue vaporization effect [9, 10, 13, 22, 51], with temperatures above 99 °C, 3) The contact surface between the blood and the tissue is smooth for flow modeling, 4) cardiac tissue is a porous medium saturated with blood, 5) deformation of porous tissue occurs due to temperature gradient, 6) local thermal equilibrium exists between the tissue and the blood flowing through it, leading to both the tissue and the blood being at one temperature at a location, 7) permeability and porosity values are isotropic and for non-infarcted tissue [22, 52], 8) blood convection occurs at the interface of blood and tissue, 9) convective heat transfer between the blood and the cardiac tissue can be described by a heat transfer coefficient, 10) initial tissue temperature is 37℃, the normal body temperature, and 11) electrical conductivity and thermal conductivity values are functions of temperature.

Electrical analysis

For radiofrequency heating near 500 kHz in RFCA, a quasistatic electric field approximation can be used since the wavelength in the tissue far exceeds tissue thickness, leading to primarily resistive heating. The rate of heat generation, \(_\) is given by

$$_=\sigma \right|}^$$

(2)

where \(\left|\mathrm\right|\) is the magnitude of the electric field (V m−1) and \(\sigma\) is the electrical conductivity of the tissue at the particular frequency (S m−1). The electric field can be calculated from the gradient of the voltage:

$$\mathrm=-\nabla \varphi$$

(3)

where \(\varphi\) is electric potential or voltage (V). The governing equations for electric potential are given as Eq. (4) [5, 22, 53].

$$\nabla \cdot (\sigma \nabla \varphi )=0$$

(4)

Boundary condition for voltage equation

Figure 4a shows the boundary conditions for electrical analysis. A constant radiofrequency power is used. The input power is set at 20W, for which the voltage used is 25 V. A voltage boundary condition is applied at the active electrode surface, while a zero flux of the electric field is imposed on all other surfaces of the catheter body (at these surfaces, catheter is electrically insulated from the blood or the cardiac tissue).

Fig. 4

Detailed boundary conditions for a voltage, b blood flow and thermal, and c poromechanical deformation

$$n\cdot \left(_-_\right)=0$$

(5)

The voltage on the bottom surface is set to 0 V to model as a dispersive electrode.

Blood flow and heat transfer analysis

Governing equations for blood velocity and temperature profiles in the cardiac tissue treats the tissue as a porous medium. The transient momentum equations (Brinkman extended Darcy model) and transient energy equations are used for analyzation [54].

Heat conduction equation

Conduction heat transfer in the catheter and electrode is given by:

$$^T}^}+\frac^T}^}+\frac^T}^}\right)}_=\frac_}\frac_}$$

(6)

where \(\alpha\) is the thermal diffusivity of the catheter or the electrode.

Momentum equations

Blood flow in the blood chamber is described by the continuity equation (Eq. 7) and the Navier–Stokes equations (Eq. (8)).

$$\frac+\frac+\frac=0$$

(7)

$$\begin\rho_b\left(\frac+u\frac+v\frac+w\frac\right)=-\frac+\mu\left(\frac+\frac+\frac\right)\\\rho_b\left(\frac+u\frac+v\frac+w\frac\right)=-\frac+\mu\left(\frac+\frac+\frac\right)-\rho_bg_y\alpha_b\left(T-T_b\right)\\\rho_b\left(\frac+u\frac+v\frac+w\frac\right)=-\frac+\mu\left(\frac+\frac+\frac\right)\end$$

(8)

Blood flow in porous tissue is described by the Brinkman extended Darcy equation (Eq. (9)) [55] that has also been used for muscles near an arteries [21]. The buoyancy term is included to account for thermally driven natural convection of blood flow:

$$\begin_\left(\frac+u\frac+v\frac+w\frac\right)=-\frac+\frac\left(\frac^u}^}+\frac^u}^}+\frac^u}^}\right)-\fracu\\ _\left(\frac+u\frac+v\frac+w\frac\right)=-\frac+\frac\left(\frac^v}^}+\frac^v}^}+\frac^v}^}\right)-\fracv-___\left(T-_\right)\\ _\left(\frac+u\frac+v\frac+w\frac\right)=-\frac+\frac\left(\frac^w}^}+\frac^w}^}+\frac^w}^}\right)-\fracw\end$$

(9)

where u, v, and w are the blood velocity components (m s−1), subscript b represents the blood phase, \(\rho\) is density (kg m−3), p is the pressure (Pa), \(\mu\) is kinematic viscosity of blood (kg m−1 s−1), \(\varepsilon\) is the tissue porosity, which is the ratio of the volume fraction of the vascular space, K is the permeability (m2), \(_\) is thermal expansion coefficient (\( }}^\)), and \(g\) is gravity (m s−2). The flow equations for the blood layer are derived from Eq. 9 as a special case of \(\varepsilon =1 \mathrm \ K=\infty\), when they revert to the Navier–Stokes equations.

The hydraulic permeability of the matrix is estimated as 3.5166 × 10–11 m2, using the Eq. (20):

$$K=\frac\sum \Delta \kappa ^$$

(10)

where K is the permeability (m2), \(\tau\) is the tortuosity (\(\tau =\sqrt)\), \(\Delta \kappa\) is the volume fraction of pores, \(\varepsilon\) is the tissue porosity, and r is the radius of pores or vessels within tissue (\(\mu\)m). The capillary diameter, d, of cardiac tissue is estimated as 30 \(\mu\)m [56]. A tissue porosity of 0.1875 is used, which is estimated using the hydraulic diameter Eq. [22, 57, 58]:

where Sv is a specific surface area (25,000 m−1) [22, 59]. In this study, tissue porosity is due to the presence of the blood vessels whose sizes (30 m\(\mu\) and 5 m\(\mu\)) are used to compute porosity and permeability, which are shown in Table 1.

Energy equations

Local thermal equilibrium can serve as a good approximation for the temperature field for certain applications involving blood vessels of small sizes [21]. The radiofrequency heat generation, as well as metabolic heat, are considered. The energy equation for the blood layer is given by:

$$_\right)}_\frac_}+_\right)}_+v\frac+w\frac\right)}_=_^T}^}+\frac^T}^}+\frac^T}^}\right)}_+_$$

(12)

while energy equations for porous cardiac tissue layer are given by:

$$_\right)}_\frac+_\right)}_\left(u\frac+v\frac+w\frac\right)=_\left(\frac^T}^}+\frac^T}^}+\frac^T}^}\right)+_+_$$

(13)

$$_\right)}_=\left(1- \varepsilon \right)_\right)}_+ \varepsilon _\right)}_\mathrm _=\left(1- \varepsilon \right)_+ \varepsilon _$$

(14)

In these equations, T is the temperature, \(_\) is the heat capacity, \(k\) is the thermal conductivity. Subscripts eff, s and b represent the effective value, solid, and blood phases, respectively, and u, v and w are the velocities in x, y, and z directions, respectively. The metabolic heat generation rate \(_\), is 684 W m−3 [60], which is the basal metabolic rate at a muscle in the thorax, that is the number of calories the body needs to accomplish its most basic life-sustaining functions. The radiofrequency heat source, \(_\), is equal to the resistive heat generated by the electric field (Eq. (2)).

Boundary condition for blood flow and heat transfer

As shown in Fig. 6b, the boundary temperatures of the porous cardiac tissue domain as well as the blood, domain are fixed at 37 °C. The outer surface between the catheter and the blood domain is considered an adiabatic boundary condition:

$$\widehat.\left(_\nabla T\right)=0$$

(15)

In the blood domain, the inlet velocity of blood in the x-direction, Uinlet, is 3 cm s−1 [61], and the outlet pressure is set to zero. For the porous tissue domain, Uinlet is assumed to be 2 m\(\upmu\) s−1 and 0.72 m\(\upmu\) s−1 corresponding to the same blood flow rate within capillary diameters of 30 m\(\upmu\) [56] and 5 m\(\mu\) [62], respectively. The outlet pressure is set to zero, as for the blood domain. At blood-catheter and blood-electrode interfaces and the remaining boundaries of the blood and porous cardiac tissue domains, no-slip conditions are applied.

Thermophysical properties

When temperature reaches 100 °C, water boils and tissue vaporizes with production of a steam pop incident [4]. The evaporation effects should be taken into consideration in this study. In medical terms, the incidence of steam pop may produce complications after the ablation procedure. Therefore, it should be prevented for effectiveness of ablation surgery outcome. The enthalpy method phase change related to temperature [9, 10, 13, 22, 51] is used for analysis, as given in Eq. (16). The density and heat capacity of both phases (liquid and gaseous) of cardiac tissue are denoted in Table 1.

$$\left\_\right\}_blood\,within\\ porous\,tissue\end} = \left\ _ _ & 0 < T \le 99^\mathrm\\ _/\Delta _ & 99 < T \le 100^\mathrm\\ __ & T>100^\circ \end\right.$$

(16)

where \(_\) and \(_\) are density and heat capacity of cardiac tissue before and post phase change, i = l as liquid phase, and i = g as gas phase. Water vaporization latent heat, \(\lambda\), is 2257 (kJ kg−1). Tissue water content, \(_\), inside cardiac tissue is 75%. Water density at 100 °C is 958 (kg m−3) and \(\Delta _\) is represented as the temperature difference as assumed by the enthalpy method. Cardiac tissue in liquid phase, \(\rho\) is 1060 (kg m−3) and \(_\) is 3111 (J kg−1 K−1), while the gas phase, \(\rho\) is 370 (kg m−3) and \(_\) is 2156 (J kg−1 K−1).

Correspondingly, the result of tissue vaporization would eventuate in desiccation. The desiccation is the dehydration effect that occurred when the cells lose water through the thermally damaged cellular wall. It results in rapid impedance (resistance) increase and then causes lower electrical conductivity. This occurrence will limit thermal volume and lead to less thermal diffusion [34]. The function of electrical conductivity, \(^,\) and thermal conductivity, \(^\), vary by temperature and can be written as temperature dependent functions as given by [9, 10, 13, 22];

$$^(T)=\left\}^, & 0 < T\le 100^\circ\\ 1.371-0.274\left(T-100\right), & < T \le 10^\circ\\ ^, & T >5^\circ \end\right.$$

(17)

$$^(T)=\left\0.531+0.0012\left(T-37\right), & 0 < T\le 100^\circ \\ 0.606, &T > 100^\circ \end\right.$$

(18)

Poromechanical deformation analysis

To obtain deformations in the cardiac tissue, a simplified quasi-static poromechanical deformation analysis is used, treating the tissue as a porous medium. The equilibrium equations for solid mechanics, written in a Cartesian coordinate system, are [63]:

$$\begin\frac_}+\frac_}+\frac_}=0\\ \frac_}+\frac_}+\frac_}=0\\ \frac_}+\frac_}+\frac_}=0\end$$

(19)

The stress–strain relationship (Eq. (20)) and the strain–displacement relationship (Eq. (21)) are as follows [63]:

$$\begin^}_=\frac\left[_-\nu \left(_+_\right)\right]+^\\ ^}_=\frac\left[_-\nu \left(_+_\right)\right]+^\\ \begin^}_=\frac\left[_-\nu \left(_+_\right)\right]+^\\ ^}_=_\left(1+\nu \right)/E\\ \begin^}_=_\left(1+\nu \right)/E\\ ^}_=_\left(1+\nu \right)/E\\ ^}_=\frac_},^}_=\frac_},^}_=\frac_}\end\end\end$$

(20)

$$\begin^}_=\frac\left(\frac_}+\frac_}\right)\\ ^}_=\frac\left(\frac_}+\frac_}\right)\\ ^}_=\frac\left(\frac_}+\frac_}\right)\end$$

(21)

where \(\xi\) denotes the stress (Pa), \(^\) is the mechanical strain, E is Young’s modulus (Pa), v is the Poisson’s ratio, and u is the average displacement (m). The thermal strain, \(^\) was calculated as follows:

where \(_\) =37 \(^\circ\) is the reference temperature and \(_\) is the temperature-dependent thermal expansion coefficient (\(^\circ\)−1) of porous cardiac tissue of temperature at 37 \(^\circ\) as indicated in Table 1. The change in porosity after deformation is defined as follows:

$$\delta =\varepsilon -_$$

(23)

$$\delta =\frac_-_^}$$

(24)

where \(\delta\) denotes as change in porosity, \(\varepsilon\) is porosity, \(_\) is volume of pore (m3), and \(d\Omega\) is total volume of solid matrix (m3). Superscript 0 represents as initial value.

Boundary condition for poromechanical deformation

As shown in Fig. 4(c), poromechanical deformation is considered only for the cardiac tissue. The boundary condition is set as moving surface for surrounding surfaces of the cardiac tissue domain. Top and bottom surfaces of tissue domain, as well as electrode-tissue interfaces, are fixed. The tissue deforms due to the thermal strain from temperature changes within the tissue domain. There is no contact force in this study. The initial stress and strain are set to zero.

Numerical computation

The electrical, heat transfer, blood flow and mechanical deformation analysis (Eqns. (2,3,4,5,6,7,8,9,10,11,12,13,14,15,16,17,18,19,20,21,22,23,24) are numerically solved using the finite element method, as implemented in the COMSOL™ Multiphysics software. The 3-D RFCA model is discretized using pyramid elements. A mesh convergence is performed to identify the suitable number of elements required as demonstrated in Fig. 5. This convergence test leads to a grid with approximately 70,000 elements.

Fig. 5

Mesh convergence finite element numerical computation showing sufficient number of elements were used in the computation