Dismiss Notice
Join Physics Forums Today!
The friendliest, high quality science and math community on the planet! Everyone who loves science is here!

BioHeat Equation solution in MATLAB using pdepe

  1. Dec 20, 2015 #1
    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:
    Code (Matlab M):

    %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
    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: Dec 21, 2015
  2. jcsd
  3. Dec 21, 2015 #2

    DrClaude

    User Avatar

    Staff: Mentor

    Why do you say it failed?

    Also, I have added [code] [/code] tags around your code to make it more readable. Please do that every time you post code.
     
  4. Dec 22, 2015 #3
    thanks for your attention.
    because It doesn't agree with my physics knowledge of this problem
     
Know someone interested in this topic? Share this thread via Reddit, Google+, Twitter, or Facebook




Similar Discussions: BioHeat Equation solution in MATLAB using pdepe
  1. MatLab pdepe (Replies: 4)

Loading...