Code covered by the BSD License  

Highlights from
High intensity focused ultrasound simulator

High intensity focused ultrasound simulator

by

 

Simulates high intensity focused ultrasound beams and heating effects in layered media

KZK_operators.m
%% Authored by Joshua Soneson 2007
function[IRK1,IRK2,CN1,CN2]=make_matrices(r,R,G,dz,dr,J,k,gamma)
%% Computes second-order diagonally implicit Runge-Kutta and 
%% Crank-Nicolson operators for integrating the kth harmonic
%% equation.

% PML vector
v = zeros(J,1);
for j=1:J
  if(r(j)>R) v(j) = 16*(r(j)-R).^2; end
end
v = exp(i*pi*v/2);
e = ones(J,1);

% build finite-difference matrix:
r(1) = 1;       % prevent divide by zero error
supdiag = (v/dr + 0.5./r)/dr;
diagonal = -2*v/dr/dr; 
subdiag = (v/dr - 0.5./r)/dr; 
A = spdiags([supdiag diagonal subdiag],-1:1,J,J);
A = conj(A');
A(1,2) = -A(1,1);
A = 0.25*A/G/k;

% split into real and imaginary parts:
ReA = real(A);
ImA = imag(A);
alpha = real(gamma)*speye(J);
beta = imag(gamma)*speye(J);

% build real-valued block system:
M = [ImA-alpha,-ReA-beta;
    ReA+beta,ImA-alpha];
I = speye(2*J);

% DIRK2 matrices:
IRK1 = I - dz*(1-1/sqrt(2))*M;
IRK2 = M; 

% CN matrices:
CN1 = I - 0.5*dz*M;
CN2 = I + 0.5*dz*M;

Contact us