Mathematical Modeling of Heat Distribution during Cryosurgery Han Liang Lim Venmathi GunasekaranIntroduction Cryosurgery, also referred to as cryoablation, is a surgical technique where undesirable or diseased tissue is frozen down using extremely low temperatures. In practice since the 1940s, Cryosurgery began as early surgeons were interested in the anesthetic properties of low temperatures. It was eventually discovered that most tissues will begin undergoing necrosis between -15°C to -40°C, and undesirable tissues were frozen from the exterior. It was not until the 1960s that the invention of surgical probes internally perfused with liquid nitrogen allowed for the insertion of such probes into the target tissue to freeze them from within. Several advantages of cryosurgery include the low invasiveness of the procedure, minimal blood flow, localizing of the site of surgery and reducing the recovery time and hospitalization time for the patient. In some instances, local anesthesia can be used in place of general anesthesia, which will result in less surgical complications. This means that in general the procedure will reduce costs for the patient. 1 While the procedure can be easily performed, it is difficult to monitor the temperature of the tissue in real time, since it would require the insertion of temperature probes which increases the invasiveness of the surgery. Thus, many researchers have turned to mathematical modeling to better understand the temperature profile of the tissue. Here, we modified the bioheat equation to obtain, where u is the temperature in the tissue, r is the radius of a sphere away from the probe, k is the thermal conductivity of tissue, C is the specific heat of tissue and Qm is the metabolic heat generation. Fig 1. A Schematic of how cryosurgery is performed on a pancreatic lesion.2Problem Formulation Establishment of geometry, boundary and initial conditions. In this study, we simulated a spherical mass of tissue that is isotropic in its thermal properties. In the middle of the tissue nests an ideal probe that occupies an infinitesimally small volume perfused by liquid nitrogen, thereby keeping the center of the sphere at -196°C. Our subsequent boundary condition is at a region infinitely far away from the probe, where the temperature should be that of body temperature, kept constant at 37°C. Before the start of the surgery, the temperature of the tissue should be constant at 37°C. However, in the freezing process, the cells will undergo a phase change at the freezing point, where they will be losing latent heat of freezing and temperature change in these cells should, theoretically, be 0. However, it has been observed clinically that the freezing state takes place across the temperature range -1°C > u > -8°C. Another point for consideration is that tissues have different thermal properties in their frozen/freezing/unfrozen states, which will be listed in the table of constants in the bioheat equation section below. Hence, our analysis will be broken down into three different temperature ranges, 37°C > u > -1°C when cells are unfrozen, -1°C > u > -8°C when cells are freezing and -8°C > u > -196°C when cells are frozen. This yields us three separate equations with different boundary conditions. Knowing the total amount of latent heat required, we took an average latent heat and combined it with the specific heat capacity constant so as to simplify the problem. We believe this range of temperature is as such since different tissues will contain a different composition of matrix components, organelles and solutes in it which will depress the freezing point. Bioheat Equation The temperature profile in the tissue can be described with Penne’s bioheat equation, which is a second ordered differential equation that goes by the form: Where C, Cb are the heat capacity of biological tissue and blood, X contains the Cartesian Coordinates x, y and z; k is the thermal conductivity of tissue, ωb is the perfusion of blood, Ta is the arterial temperature, u is the tissue temperature and Qm is the metabolic heat generation. The values are presented in the table of constants below.This equation can be further simplified in our instance if we consider that in rapidly freezing tissues, we will first cause vasoconstriction in the capillaries before freezing all the blood in the capillaries. In the absence of perfusion, the already small ωb term goes to zero. Also, cells will not be able to generate any metabolic heat when frozen, and Qm is nonexistent in temperatures below 0. Putting these together, before the cells are freezing we have a nonhomogenous differential equation and cells in the frozen and freezing states can be described with a homogenous differential equation instead. Upon inspection, it becomes apparent that the diffusion in the system is radial in spherical coordinates and independent of the other spherical coordinates, φ and θ. By converting Penne’s bioheat equation into spherical coordinates, we obtain the equation first mentioned in our introduction, Which we will solve analytically to obtain our temperature profile. It is worth mentioning here that in converting from a Cartesian coordinate system, we generate an additional “ ” term describing the dependence of change in temperature with time on the spatial variation of temperature in both the first and second order, which hints at a solution utilizing spherical Bessel’s functions. Analytical Solution Starting with our modified bioheat equation, 22()mu u uc k Qt r r r         …(1) Fig 2. A table of constants for values used in the bio heat equation for all three temperature rangesAnd our boundary and initial conditions, 1(0, )u t T and ( , )bu b t T. (2 and 3) 1181728ufuuc u CccQc C u Cc u C      (4) , 18128ufuuk u Ckkk C u Ck u C     (5) 0.042 10 8 108muCQ C u CuC     (6) 137 1( ,0) 1 8 188uCu r T C u CuC         (7) Where (2) and (3) are our boundary conditions, (7) is our initial condition and (4) through (6) are our constants. To solve the non-homogeneous equation