Matlab based temperature solver

From MGH/MIT Parallel Transmission Resources
Jump to navigation Jump to search

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.