BioHeat Equation solution in MATLAB using pdepe

  • Context: MATLAB 
  • Thread starter Thread starter m.r.fouladi
  • Start date Start date
  • Tags Tags
    Matlab Pde
Click For Summary
SUMMARY

The forum discussion centers on solving the bioheat equation using MATLAB's pdepe function. The equation is defined as ∂T/∂t = α ∇2T + 1/ρc[S+Sp+Sm], with specific constants provided for a tumor located at the center of healthy tissue. The user seeks assistance in correcting their MATLAB code, which is intended to model temperature variation in the tumor due to energy delivered by an electrode. The provided code has issues that prevent it from aligning with the physical principles of the problem.

PREREQUISITES
  • Understanding of the bioheat equation and its parameters
  • Familiarity with MATLAB and the pdepe function
  • Knowledge of numerical methods for partial differential equations
  • Basic principles of heat transfer in biological tissues
NEXT STEPS
  • Review MATLAB documentation on the pdepe function for solving PDEs
  • Study the implementation of boundary conditions in MATLAB PDEs
  • Learn about numerical stability and convergence in finite difference methods
  • Explore the physical implications of the bioheat equation in medical applications
USEFUL FOR

Researchers, biomedical engineers, and students involved in computational modeling of thermal effects in biological tissues, particularly those using MATLAB for simulations.

m.r.fouladi
Messages
2
Reaction score
0
We have this Equation as bioheat equation:
∂T/∂t = α ∇2T + 1/ρc[S+Sp+Sm]
and also this:
Sp=mbcb(Tab-T)
that all α,ρ,c,S,Sm,mb,cb,Tab are constants, now I want to solve this equation in conditions below with pdepe in MATLAB:
There is a Tumor as a sphere with radius 1 cm exactly in center of a Normal Tissue with radius of 5 cm, an electrode at t=0 gives an Energy to the center of Tumor for 400 seconds that T will vary and we want to have this variation by solving this equation, also below assumptions is considered:
We assume that Tumor has the radius of 1 cm and that is perfused at the same level as healthy tissue. The tumor is located at the center of a sphere of healthy tissue that has a radius of 5 cm and all of the tissue has a metabolic energy release rate of 145 W/m3.
The power is delivered locally to the tumor by positioning the electrode in the center of the tumor. the initial temperature of tumor and surrounding tissue is uniform at 37 °C. in addition, It is assumed that because of symmetry there is zero flux at the center of the tumor and that at the outer surface of healthy tissue the flux is also zero. the following numerical values for the parameters are used:
mb = 0.18
cb = 3300
Sm = 145
Tab = 37 °C
α = 10∧(-7)
c=3800
S= 4 * 10∧(5)
ρ= 850.I coded it in MATLAB such below, but my code isn't right as my knowledge, please correct it my masters:
Matlab:
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
function pdex1

m = 2;
x = 0.001:0.001:05;
t = 0:1:400;

sol = pdepe(m,@pdex1pde,@pdex1ic,@pdex1bc,x,t);
% Extract the first solution component as u.
u = sol(:,:,1);

% A surface plot is often a good way to study a solution.
surf(x,t,u)
title('Numerical solution computed with 20 mesh points.')
xlabel('Distance x')
ylabel('Time t')

% A solution profile can also be illuminating.
figure
plot(x,u(end,:))
title('Solution at t = 2')
xlabel('Distance x')
ylabel('u(x,400)')
% --------------------------------------------------------------
function [c,f,s] = pdex1pde(x,t,u,DuDx)
c = 1e7;
f = DuDx;
if (x<0.01)
s = (18.089*10^(-2)-(18*10^(-5)*u))*(1e7);
else
s=594*(310-u)/(850*3800)
end
% --------------------------------------------------------------
function u0 = pdex1ic(x)
u0 = 310;
% --------------------------------------------------------------
function [pl,ql,pr,qr] = pdex1bc(xl,ul,xr,ur,t)
pl = ul-exp((-18)*10^(-5)*t)-1004.94;
ql = 0;
pr = ul-exp((-18)*10^(-5)*t)-1004.94;
qr = 0;

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
why It's failed?
 
Last edited by a moderator:
Physics news on Phys.org
m.r.fouladi said:
why It's failed?
Why do you say it failed?

Also, I have added
Code:
tags around your code to make it more readable. Please do that every time you post code.
 
DrClaude said:
Why do you say it failed?

Also, I have added
Code:
tags around your code to make it more readable. Please do that every time you post code.
thanks for your attention.
because It doesn't agree with my physics knowledge of this problem
 

Similar threads

  • · Replies 8 ·
Replies
8
Views
3K
Replies
4
Views
2K
  • · Replies 6 ·
Replies
6
Views
4K
  • · Replies 1 ·
Replies
1
Views
5K
  • · Replies 1 ·
Replies
1
Views
4K
  • · Replies 9 ·
Replies
9
Views
3K
  • · Replies 3 ·
Replies
3
Views
4K
  • · Replies 7 ·
Replies
7
Views
4K
  • · Replies 6 ·
Replies
6
Views
2K
  • · Replies 2 ·
Replies
2
Views
5K