Modelling Of Heat Conduction Along Vein Walls In Endovenous Laser Treatment
W Fuller, M Raiti, R Bush
Citation
W Fuller, M Raiti, R Bush. Modelling Of Heat Conduction Along Vein Walls In Endovenous Laser Treatment. The Internet Journal of Bioengineering. 2010 Volume 5 Number 1.
Abstract
Laser ablation of the saphenous vein uses laser-tipped probes to produce photothermal effects in the vein. In this study we construct a mathematical model of the effects of laser-induced thermal heating and conduction along the vein wall. We formulate and solve the relevant two-dimensional heat conduction problem. Simulations of the solution to the problem resolve an aspect of a debated question of the procedure.
Introduction
Since its inception, controversy has surrounded the exact etiology of thermal damage to the saphenous vein [1-5]. There are proponents of both conduction energy (laser tip contact) and convective energy (steam bubbles) as being the etiology of resultant tissue damage. This controversy has been debated not only from the scientific aspect, but in the legal arena as well.
As reported in previous publications [1, 5], the hallmark of acute damage at the time of thermal ablation is complete loss of endothelium. This finding is present when laser firing at 12-15 watts is done at intervals of 3-4 mm apart. If it is assumed that the laser point contacts the surface of the vein, which occurs when a conventional bare tip fiber is used, then the question arises of whether the contact and the resultant conduction by tissue is responsible for the histologic findings. The absence of endothelium immediately following laser thermal ablation is a constant finding regardless of wavelength and pulse duration. Based on our histologic observations it is our hypothesis that heat conduction would not be the primary etiology of these observed luminal changes. The purpose of this article is to present the results of testing our hypothesis by simulations based on the construction of a relevant mathematical model which will accurately calculate both temperature and distance traveled of energy from the laser application site.
Construction of the model
The key parameter of the model is the initial temperature to which the laser probe elevates the affected region of the vein wall. In the model we will keep this temperature as a free parameter, but when we present results of simulations based on the model, we will assume that at the time of laser firing the initial temperature at the laser tip is 100 C. This assumption is reasonable since (1) steam bubbles occur at this point, and (2) we have measured the temperature of these bubbles impinging on the vein wall. We have run simulations at other initial temperatures, such as 700 C and the associated results do not significantly change our final conclusions.
The first step in building a valid mathematical model is the enunciation of the assumptions on which the model rests. We base our model on five postulates derived from both experience and theory.
Discussion: The rectangular shape is chosen so as to model the region where the probe tip contacts the partially collapsed vein wall. Our analysis will indicate that the time evolution of the heat conduction smooths out the boundary of the initial region so that the assumption of a rectangular shape is not of particular importance in determining the region of cell death. As indicated above, we will present the results of simulations for a value of T0 of 100 C.
Discussion: That the saphenous vein is a uniform right circular cylinder is a simplifying assumption that models the actual vein as it appears in ultrasounds of the procedure and in clinical experience. The tumescent solution outside of the vein certainly exerts no tensile forces on the vein wall, and any compression it might provide has minimal effect on the surface area. The assumption of infinite length merely indicates that the regions of laser application are far enough from the ends of the vein that the longitudinal boundaries play no role in the local heat conduction process.
Discussion: This postulate, along with Postulate 2, restricts the phenomenon at hand to one of heat conduction along the vein wall. The assumption that heat conduction is the primary mechanism causing endothelial necrosis thus requires Postulates 2 and 4 as its underlying support. From this point of view the effect of energy transfer through the blood and into the vein wall is a separate problem. It is true that dissipation of energy into the ambient anesthetic and tissue occurs. Experimental results [6] indicate that such heat loss is minimal. In any case assumption of Postulate 4 ensures that whatever region of cell death we obtain in this model will be an upper bound for the actual region. Any model incorporating other modes of energy loss will certainly produce smaller regions of cell death.
T = T(x, y, z, t) is the temperature at the point with coordinates
Discussion: Postulates 4 and 5 imply that heat flows along the vein wall in such a way that energy is conserved. Standard energy conservation arguments indicate that the local behavior of the heat flow is described by the heat conduction equation [7].
Here ΔT is the Laplacian of T and the constant k is the thermal diffusivity of the vein wall. We take for the value of thermal diffusivity that of human myocardium tissue: k = 1.289 x 10-7 m2/sec.
The heat conduction problem: the planar approximation
Consider a right-handed coordinate system with the z-axis along the central longitudinal axis of the vein and the x-axis through the center of the rectangle initially raised to a temperature T0 on the vein wall. In cylindrical coordinates equation (1) becomes
Here φ is the azimuthal angle measured counterclockwise from the x-axis and ρ is the radial distance out from the origin. Along the vein wall ρ is a constant, R, so that on the cylindrical surface the equation becomes
We take the value of the radius of the vein as ρ = R = 4 mm.
We may shift to a local coordinate system on the cylinder and centered on the rectangle initially at temperature T0. The local coordinates then become the circumferential distance
s = Rφ and the longitudinal distance z as shown in Figure 1. Also shown is the initial condition rectangle in red. The length 2a is the diameter of the probe tip, while 2b is the length of the probe tip. Typically a = .25 mm and b = 2 mm.
In these local coordinates equation (3) assumes the form:
The initial condition problem associated with equation (4) has circumferential boundary conditions which are periodic with period 2pR. Since equation (4) is a diffusion-type equation, we expect the temperature distribution to spread out from its initial rectangular region. When it spreads as far as ±pR ~ 12.57 mm there will therefore be an overlapping of the left and right traveling “tails” of the solution. The linearity of equation (4) ensures that the overlap solution is the superposition of the tails. On the other hand, for the time periods of importance overlap effects are very small and, as we shall see, the spreading heat distribution is essentially contained within the region - pR
Here Θ(x) is the Heaviside function.
Solution of the planar approximation problem
Our method of solution will be to construct the two-dimensional Green’s function for equation (4) and then to integrate the product of the Green’s function and the initial condition.
[See for example 8, 9]
The relevant Green's function problem is to solve equation (4) subject to the initial condition, T(
(
Taking the Fourier transform of (4) yields an ordinary differential equation:
Here a2 =
Solving for τ gives
Taking the inverse transform yields the Green's function:
The initial condition problem has solution:
Upon substitution of the initial condition this becomes:
The solution we have obtained is the product of two independent Gaussian distributions spreading along the vein wall. One Gaussian spreads longitudinally while the other spreads azimuthally. The formula obtained for T(s, z, t) permits the simulation of the conduction of heat away from the laser tip along the vein wall. Figures 2 and 3 show the results of such a simulation for an initial temperature of T0 = 100 C and at times t = .63 sec and t = 5.00 sec. In each figure the region of possible cell death is indicated in white. We estimate that, in each 15 mm longitudinal segment of the saphenous vein, not more than 14 per cent of the vein wall experiences temperatures above the cytotoxic threshhold. In the figures one can also observe the smoothing out of the boundary of the initial condition rectangle. Such smoothing of boundaries under heat conduction is well-known and indicates that in our Postulate 1 other assumptions on the shape of the initial temperature distribution would give results similar to ours.
Figure 14
Figure 15
Figure 4 displays the azimuthal Gaussian distribution for the temperature ratio, R = T/T0, shortly after the burst of energy from the laser probe at t = 1.25 sec. Figure 5 shows the same distribution at a later time t = 10.0 sec. Since overlapping of the tails of the temperature distribution occurs for s = ± .000785 m, these graphs indicate that for small times there is no appreciable overlap effect. At large times, however, the overlap effect may be as much as twenty per cent of the maximum temperature at the corresponding time. In the next section we shall construct a solution that will be valid at these longer times.
Figure 16
The initial condition heat conduction problem with periodic boundary condition
The problem at hand is to incorporate the circumferential boundary conditions into the initial condition problem; consequently, we seek a solution which is periodic in s with period 2πR. In this case the problem may be formulated as:
The sequence of functions
for n belonging to , forms an orthonormal set on the interval [-πR, + πR]. We expand T(t, s, z) in terms of this basis:
This Fourier series expansion for T(t, s, z) manifestly satisfies the periodic boundary condition.
Substitution of the expansion (7) for T(t, s, z) into the partial differential equation (6) yields
The linear independence of the functions {en (s)} gives for each n:
The initial condition dependent on z is specialized from f(s, z) and may be written:
where is a constant to be determined. Let σn (t, ω) be the Fourier transform of τn (t, z) with respect to z. The Fourier transform of equation (8) yields an ordinary differential equation for σn which has solution:
Here σ0n is a constant determined by the initial condition for τn (t, z).
The Green’s function solution for τn (t, z) may be obtained by assuming that for each n τn (0,z) = δ(z – ζ). σn (t, ω) accordingly becomes
Taking the inverse Fourier transform yields the Green’s function solution:
The solution for τn (t, z) may be obtained by integrating Гn (z|ζ, t) against the initial condition (9). Substituting the resulting expression for τn (t, z) back into equation (7) produces the general solution for T(t, s, z) with t > 0:
The coefficients τ0n may be determined by use of the azimuthal initial condition which requires
The orthogonality of the basis functions {en (s)} permits calculation of the τ0n :
and for n ≠ 0:
Substitution into equation (10) provides the desired solution for the temperature function:
{image:31}
Figures 6 and 7 show the results of simulations of this solution for an initial temperature of T0 = 100 C and at times t = .63 sec and t = 5.00 sec. In both simulations the series in s are expanded to include the first 1500 terms. As before the region of possible cell death is indicated in white. The regions above the cytotoxic limit are smaller than the regions at the same times in the simulation results given in the previous section. The similarity of the results thus indicate the validity of the planar approximation approach.
{image:32}
{image:33}
Discussion
Ultrasounds taken during endovenous ablation consistently show the formation of steam bubbles as demonstrated by Proebstle et al. [10] and Yilmaz et al. [3]. Since the erythrocyte is 80% water and considering the affinity of these wavelengths for water, it is reasonable to assume that these bubbles represent the formation of steam.
There are three types of energy that can account for the observed tissue damage. These are conductive (heat transfer into the target tissue), convective (heat transfer into the liquid), and radiant (light absorption). Much of the radiant energy is rapidly absorbed by the heme protein and water. The model we have examined in this paper is based on the assumption that all of the energy from the laser probe results in conductive energy occurring at the point of laser contact. This assumption has provided a way of testing the most extreme case of conductive damage to the endothelium. In view of the results of our simulations it is reasonable to conclude that large-scale endothelial necrosis is a result of energy transfer in the convective and radiant channels.
Based on our mathematical model, tissue conduction by laser application only accounts for about 15% of cellular damage when temperatures of 45 C are used as the initial point of irreversible cell damage. This point is well below the 50 C needed for cell death. When 50 C is considered, the area of destruction would be even smaller.
Conclusion
We conclude from the preceding discussion that massive cell death of the endothelium cannot be due solely to conduction along the inner surface of the vein wall. Laser tip contact is not necessary for thermal ablation to be effective. Experimental evidence [5] indicates that the fluid in the vein flows in a highly turbulent manner and that portions of the fluid undergo a liquid-gas phase transition. The ordinary diffusion equation is no longer valid in such a regime. Modelling heat transfer in this regime remains an open question.