**3. Numerical modeling**

This section discusses the methodology for solving the RTE (Eq. 7) using the DOM. Then, the numerical solution of the DPL model-based bio-heat transfer is obtained using the FVM to determine the temperature distribution inside the laser-irradiated biological tissue.

In DOM, the RTE, an integrodifferential equation, has been transformed into a set of partial differential equations [10]. So, the RTE for the diffuse component of intensity (Eq. (12)) with the corresponding boundary conditions (Eq. 13) in any given discrete direction (*Ω<sup>m</sup>*,*<sup>n</sup>*) can be written as follows

$$\begin{split} \frac{1}{c} \frac{\partial I\_d^{m,n}}{\partial t} + \mu^{m,n} \frac{\partial I\_d^{m,n}}{\partial \mathbf{x}} + \eta^{m,n} \frac{\partial I\_d^{m,n}}{\partial \boldsymbol{\eta}} &= -\beta I\_d^{m,n} + \frac{\sigma\_s}{4\pi} \sum\_{m'=1}^{N\_\theta} \sum\_{n'=1}^{N\_\theta} I\_d^{m',n'} \mathcal{O}^{m',n';m,n} \alpha\_\theta^{m'} \alpha\_\boldsymbol{\eta}^{n'} \\ &+ \frac{\sigma\_s}{4\pi} I\_c^{m\_c,n\_c} \mathcal{O}^{m\_c,n\_c;m,n} \end{split} \tag{12}$$

*Modeling of Laser-Irradiated Biological Tissue DOI: http://dx.doi.org/10.5772/intechopen.106794*

$$\text{Left wall } (\varkappa = 0, \ 0 < \jmath < W) : I\_d^{m, n} = 0, \mu^{m, n} > 0 \tag{13}$$

$$\text{Right wall } (\varkappa = L, \ 0 < \jmath < W) : I\_d^{m, n} = 0, \mu^{m, n} < 0 \tag{14}$$

$$\text{Bottom wall } (y = 0, \ 0 < \infty < L): I\_d^{m, n} = 0, \eta^{m, n} > 0 \tag{15}$$

$$\text{Top wall } (\text{y} = \text{W, } 0 < \text{x} < L): I\_d^{m, n} = 0, \eta^{m, n} < 0 \tag{16}$$

The boundary conditions Eqs. (13)–(16) have been obtained using Eq. (8) under the assumption that the boundary is non-reflecting and the magnitude of the light intensity emitted from the surface is negligible compared to the intensity of the shortpulse laser used.

Using a fully implicit backward differencing scheme in time, Eq. (12), after simplification, becomes

$$\mu^{m,n}\frac{\partial I\_d^{m,n}(t^\*)}{\partial x} + \eta^{m,n}\frac{\partial I\_d^{m,n}(t^\*)}{\partial y} + \left[\frac{\beta}{\Delta t^\*} + \beta\right]I\_d^{m,n}(t^\*) = \frac{\beta}{\Delta t^\*}I\_d^{m,n}(t^\* - \Delta t^\*) + S\_t^{m,n} \tag{17}$$

where *<sup>t</sup>* <sup>∗</sup> <sup>¼</sup> *<sup>β</sup>ct*,

$$\mathbf{S}\_{t}^{m,n} = \frac{\sigma\_{\varsigma}}{4\pi} \sum\_{m'=1}^{N\_{\theta}} \sum\_{n'=1}^{N\_{\phi}} I\_{d}^{m',n'} \mathcal{Q}^{m',n';m,n} \alpha\_{\theta}^{m'} \alpha\_{\phi}^{n'} + \frac{\sigma\_{\varsigma}}{4\pi} I\_{c}^{m\_{\varsigma},n\_{c}} \mathcal{Q}^{m\_{\varsigma},n\_{c};m,n} \tag{18}$$

Integrating Eq. (17) over the control volume Δ*V* leads to

$$\begin{split} & \mu^{m,n} \Delta \eta \left[ I\_{d,E}^{m,n}(t^\*) - I\_{d,W}^{m,n}(t^\*) \right] + \eta^{m,n} \Delta \mathbf{x} \left[ I\_{d,N}^{m,n}(t^\*) - I\_{d,S}^{m,n}(t^\*) \right] \\ &= -\left[ \frac{\beta}{\Delta t^\*} + \beta \right] \Delta V I\_{d,P}^{m,n}(t^\*) + \Delta V S\_{t,P}^{m,n} + \frac{\beta \Delta V}{\Delta t^\*} I\_{d,P}^{m,n}(t^\* - \Delta t^\*) \end{split} \tag{19}$$

where, *I m <sup>d</sup>*,*<sup>P</sup>* is the intensity at the cell center *P*, *Sm <sup>t</sup>*,*<sup>P</sup>* is the source term at the cell center *P* and Δ*V* ¼ Δ*x* � Δ*y*. For relating the cell-surface intensities with cell-center intensity, the step differencing scheme has been used to avoid any possibility of providing unphysical results [10].

Eq. (19) is simplified using the step differencing scheme, and the resulting equation becomes as:

$$I\_{d,P}^{m,n}(t^\*) = \frac{\left[\mu^{m,n}\Delta\eta I\_{d,W}^{m,n}(t^\*) + \eta^{m,n}\Delta\eta I\_{d,S}^{m,n}(t^\*) + \Delta V S\_{t,P}^{m,n} + \left(\frac{\beta\Delta V}{\Delta t^\*}\right)I\_{d,P}^{m,n}(t^\* - \Delta t^\*)\right]}{\mu^{m,n}\Delta\eta + \eta^{m,n}\Delta\eta + \left[\frac{\beta}{\Delta t^\*} + \beta\right]\Delta V} \tag{20}$$

As earlier mentioned, the solution of RTE provides the intensity distribution, which is used to determine the divergence of the radiative heat flux (Eq. 5). Subsequently, find the temperature rise due to a single pulse by solving the Eq. (9). Then, Eq. (10) has been solved using the FVM to determine the temperature distribution inside the biological tissue subjected to the train of pulses. The algorithm for solving this multi-time scale problem has been discussed in Section 2.

Under the Cartesian coordinates system, the governing equation (Eq. 11) in *x*- directions can be expressed as [21]

$$\frac{\pi\_q}{a}\frac{\partial^2 q\_x}{\partial t^2} + \frac{1}{a}\frac{\partial q\_x}{\partial t} = \frac{\partial}{\partial \mathbf{x}} \left(\frac{\partial q\_x}{\partial \mathbf{x}} + \frac{\partial q\_y}{\partial \mathbf{y}}\right) + \tau\_T \frac{\partial}{\partial t} \left\{\frac{\partial}{\partial \mathbf{x}} \left(\frac{\partial q\_x}{\partial \mathbf{x}} + \frac{\partial q\_y}{\partial \mathbf{y}}\right)\right\} \tag{21}$$

Integrating Eq. (21) over a given control volume and over the time step of Δ*t,* i.e., from *t* to *t*+Δ*t* followed by discretization based on the backward difference in time leads to the following discretized form of the Eq. (22) [21].

$$a\_{\rm xP}q\_{\rm xP}^{\rm t+\Delta t} = a\_W q\_{\rm xW}^{\rm t+\Delta t} + a\_E q\_{\rm xE}^{\rm t+\Delta t} + b\_{\rm x} \tag{22}$$

$$\text{where } \mathfrak{a}\_{W} = \frac{\Delta y}{(\delta x)\_{w}} (\Delta t + \mathfrak{r}\_{T}), \mathfrak{a}\_{E} = \frac{\Delta y}{(\delta x)\_{e}} (\Delta t + \mathfrak{r}\_{T}), \mathfrak{a}\_{\mathbf{x}P} = \frac{\Delta x \Delta y}{a} \left(\mathfrak{1} + \frac{\mathfrak{r}\_{q}}{\Delta t}\right) + \mathfrak{a}\_{W} + \mathfrak{a}\_{E} \text{ (}$$

*bx* <sup>¼</sup> <sup>Δ</sup>*x*Δ*<sup>y</sup> <sup>α</sup>* <sup>1</sup> <sup>þ</sup> 2*τ<sup>q</sup>* Δ*t* þ Δ*yτ<sup>T</sup>* ð Þ *δx <sup>w</sup>* þ Δ*yτ<sup>T</sup>* ð Þ *δx <sup>e</sup> qx t <sup>P</sup>* � <sup>Δ</sup>*yτ<sup>T</sup>* ð Þ *δx <sup>w</sup> qx t <sup>W</sup>* � <sup>Δ</sup>*yτ<sup>T</sup>* ð Þ *δx <sup>e</sup> qx t E* � <sup>Δ</sup>*x*Δ*yτ<sup>q</sup> <sup>α</sup>*Δ*<sup>t</sup> qx t*�Δ*t <sup>P</sup>* þ ð Þ Δ*t* þ *τ<sup>T</sup> qy t*þΔ*t ne* � *qy t*þΔ*t se* � *qy t*þΔ*t nw* <sup>þ</sup> *qy t*þΔ*t sw* � *τ<sup>T</sup> qy t ne* � *qy t se* � *qy t nw* <sup>þ</sup> *qy t sw*

Here, the symbols *t+*Δ*t*, *t,* and *t-*Δ*t* represent the current time step, previous time step, and previous to previous time step, respectively.

Similarly, the discretized form of the equation for the heat flux component in the *y*-direction of Eq. (11) can be obtained.

Once the heat flux component in both directions has been obtained, the Eq. (10) has been discretized using the FVM to determine the temperature distribution at the current time instant using Eq. (23) [21].

$$T\_P^{t+\Delta t} = T\_P^t + \frac{\Delta t}{\rho c\_v \Delta x} \left( q\_{xw}^{t+\Delta t} - q\_{x\_\varepsilon}^{t+\Delta t} \right) + \frac{\Delta t}{\rho c\_v \Delta y} \left( q\_{y\_\varepsilon}^{t+\Delta t} - q\_{y\_n}^{t+\Delta t} \right) \tag{23}$$

The obtained numerical results have been compared with the analytical solution. The analytical solution of the DPL model-based bio-heat transfer equation has been obtained using the FIT technique, which is discussed in Section 4.
