# Matlab based temperature solver

### Citation

Bastien Guérin, Jorge F. Villena, Athanasios G. Polimeridis, Elfar Adalsteinsson, Luca Daniel, Jacob K. White, Bruce R. Rosen and Lawrence L. Wald (2016). "Ultimate hyperthermia: Computation of the best achievable radio-frequency hyperthermia treatments in non-uniform body models". Annual meeting of the International Society for Magnetic Resonance in Medicine. Singapore.

### What this code does

This is a Crank-Nicholson Matlab-based temperature solver with heat generation (=SAR) and tissue perfusion. The code mixed spatial and temporal finite differences using the Crank-Nicholson scheme. Non-homogeneous tissues properties (including tissue perfusion, thermal conductivity and heat capacity) are supported. The code is fully validated by comparison with the SEMCAD temperature solver. Two solve options are available: Transient and steady-state calculations.

### Download files

The code can be downloaded here.

### Instructions

#### Transient simulations

To run a simple example, go to the sub-folder SEMCAD_validation/Simu_1_Results and run the script "semcad_vs_pbhs.m". This runs a comparison of the transient solver with SEMCAD. This script should make it clear how the code is run. The input to the code are, in this order:

- T0: The initial temperature map (can be non-uniform) in K
- DUR: The total duration of the simulation in s
- [dx dy dz]: A 1x3 Matlab vector with the x, y and z spatial resolutions in m
- RHO: The density map (can be non-uniform) in kg/m^3
- C: The heat capacity map (can be non-uniform) in J/kg/K
- K: The thermal conductivity map (can be non-uniform) in W/m/K
- Wb: The blood perfusion map (can be non-uniform) in kg/s/m^3
- Cb: The blood heat capacity in J/kg/K
- Tb: The blood temperature in K
- Q: The heat rate generation map (i.e., SAR*RHO which can be non-uniform) in W/m^3

#### Steady state simulations

To run a simple example of the steady-state simulation, go to the TEST_SteadyState/ folder and run "TEST_SS.m". Again, the way to run the code should be straightforward. The input to the steady-state function are similar to the input to the transient solve function, except obviously the duration is not needed.