Pennes' bioheat equation, a cornerstone of thermal modeling in biological tissues, describes how temperature varies within living tissue. It extends the standard heat diffusion equation by incorporating the effects of blood perfusion and metabolic heat generation. The equation accounts for tissue properties like density, specific heat capacity, and thermal conductivity, along with blood perfusion rate, arterial temperature, and metabolic rate.

Weak Formulation

To apply FEM, we first derive the weak form of the equation by multiplying both sides with a test function $v$ and integrating over the domain $\Omega$ :

$$ \int_{\Omega} \rho_t c_t \frac{\partial \theta}{\partial t} v d \Omega=\int_{\Omega} \nabla \cdot(K \nabla \theta) v d \Omega+\int_{\Omega} \rho_b c_b \omega_b\left(\theta_a-\theta\right) v d \Omega+\int_{\Omega} M_A v d \Omega $$

Applying the divergence theorem to the conduction term:

$$ \int_{\Omega} \nabla \cdot(K \nabla \theta) v d \Omega=-\int_{\Omega} K \nabla \theta \cdot \nabla v d \Omega+\int_{\Gamma} K \nabla \theta \cdot n v d \Gamma $$

where $\Gamma$ represents the domain boundary and $n$ is the outward normal. Rewriting, we obtain the weak form:

$$ \int_{\Omega} \rho_t c_t \frac{\partial \theta}{\partial t} v d \Omega+\int_{\Omega} K \nabla \theta \cdot \nabla v d \Omega+\int_{\Omega} \rho_b c_b \omega_b \theta v d \Omega=\int_{\Omega} \rho_b c_b \omega_b \theta_a v d \Omega+\int_{\Omega} M_A v d \Omega+\int_{\Gamma} q v d \Gamma $$

where $q=K \nabla \theta \cdot n$ is the heat flux on the boundary.

Finite Element Discretization

We approximate the temperature field $\theta(x, t)$ using basis functions $\phi_i(x)$ :

$$ \theta(x, t) \approx \sum_{i=1}^N \theta_i(t) \phi_i(x) $$

Similarly, we choose test functions $v$ from the same space:

$$ v=\phi_j $$

Substituting these approximations into the weak form results in the discrete FEM formulation:

$$ M \frac{d \theta }{d t}+ K \theta + W \theta = F $$

where: