Discover MakerZone

MATLAB and Simulink resources for Arduino, LEGO, and Raspberry Pi

Learn more

Discover what MATLAB® can do for your career.

Opportunities for recent engineering grads.

Apply Today

New to MATLAB?

Solving one-dimensional PDE's using the PDE Toolbox

Asked by Abed Alnaif

Abed Alnaif (view profile)

on 8 Dec 2012
Latest activity Commented on by Abed Alnaif

Abed Alnaif (view profile)

on 20 Mar 2015 at 14:51

Hi,

I've been trying to solve a non-linear, heat-equation-type system of PDE's using the 'pdepe' function, with only one dimension in space. However, for many sets of parameter values, the solver exhibits unstable behaviour (oscillations, etc). I think that my problem demands a more sophisticated solver, due to nonlinearities and discontinuities. Hence, I am now trying to use the PDE Toolbox, hoping that it would be able to handle my problem, since it has an adaptive mesh algorithm, etc. However, I was unable to figure out how to use the PDE Toolbox for one-dimensional problems, nor was I able to find any examples of this.

Is it possible to use the PDE Toolbox to solve one-dimensional problems? I'm find with using the command-line interface (I don't necessarily need to use the GUI). If this is possible, I'd greatly appreciate it if someone could provide me with some code that solves a simple heat equation PDE (one dimension in space) using the PDE Toolbox, just to get me on the right foot (since, for some reason, I get really lost when I try to read the documentation for the PDE Toolbox).

Thanks,

Abed

0 Comments

Abed Alnaif

Abed Alnaif (view profile)

Products

No products are associated with this question.

2 Answers

Answer by Bill Greene

Bill Greene (view profile)

on 8 Dec 2012
Accepted answer

Hi,

Solving a 1-D PDE with PDE Toolbox is fairly straightforward. You just define a rectangular region of the appropriate width and arbitrary height. On the two edges at y=constant, you want a zero-Neumann BC. On the edges at x=constant, the BCs should not vary in y.

I've appended a very simple example of time-dependent heat transfer in a bar below.

Based on your description of the problem you are trying to solve, however, I can't think of any reason why PDE Toolbox should give a better solution than pdepe.

Bill

function simple1DTest
% 1D transient heat transfer in x-direction
h =.1; w=1; % width equal one, height is arbitrary (set to .1)
g = decsg([3 4 0 w w 0 0 0 h h]', 'R1', ('R1')');
[p, e, t]=initmesh(g, 'Hmax', .05);
b=@boundFile;
c=1; a=0; d=1;
f='100*x'; % heat load varies along the bar
u0 = 0; % initial temperature equals zero
tlist = 0:.02:1;
u=parabolic(u0, tlist, b,p,e,t,c,a,f,d);
figure; pdeplot(p,e,t, 'xydata', u(:,end), 'contour', 'on'); axis equal;
title(sprintf('Temperature Distribution at time=%g seconds', tlist(end)));
figure; plot(tlist, u(2,:)); xlabel 'Time'; ylabel 'Temperature'; grid;
title 'Temperature at tip as a function of time'
end
function [ q, g, h, r ] = boundFile( p, e, u, time )
N = 1; ne = size(e,2);
q = zeros(N^2, ne); g = zeros(N, ne);
h = zeros(N^2, 2*ne); r = zeros(N, 2*ne);
% zero Neumann BCs (insulated) on edges 1 and 3 at y=0 and y=h
% and the right edge, edge 2
for i=1:ne
  switch(e(5,i))
    case 4
      h(1,i) = 1; h(1,i+ne) = 1;
      r(i) = 500; r(i+ne) = 500; % 500 on left edge
  end
end
end

1 Comment

Abed Alnaif

Abed Alnaif (view profile)

on 10 Dec 2012

Thanks, I'll give it a shot.

Bill Greene

Bill Greene (view profile)

Answer by tom_brot

tom_brot (view profile)

on 20 Mar 2015 at 11:16
Edited by tom_brot

tom_brot (view profile)

on 20 Mar 2015 at 11:18

Ok this question has already been answered, but I want to share my experience with pdepe.

I also had problems with oscilating, unstable behaviour of pdepe. I found that the Problem lies not directly in the pdepe, but in the underlying ode15s solver. You can specify some options (with odeset()) for the pdepe-solver which are then passed to the ode15s solver.

The problem is, that only very few of the possible odeset() options are passed. To get a much more stable solution from the underlying ode15s solver I reduced the 'MaxOrder' option from the Default value of {5} to 2.

It's not that hard to make the option available for your pdepe solver.

How to make your pdepe solver more stable:

1) Find the file 'pdepe.m' in somewhere in your matlab installation path and copy it to your working directory. And rename it to e.g. 'pdepe_opt.m'.

2) Find the file 'pdentrp.m' (which is located in the 'private' directory next to 'pdepe.m') and copy it to your working directory.

3) In your local copy of pdepe.m, namely pdepe_opt.m add the line

maxorder = odeget(options,'MaxOrder',[],'fast');

after the line

maxstep = odeget(options,'MaxStep',[],'fast');

4) Then add

'MaxOrder',maxorder,

as argument of the odeset command, which is some lines below.

5) Before the line in your code where you call pdepe, insert

options = odeset('MaxOrder',2)

(or any number from 1 to 5) and then call the optimized pdepe function with your specified options:

sol = pdepe_opt(m,@pdex1pde,@pdex1ic,@pdex1bc,x,t,options)

I hope this works as well for you as it did for me.

Best regards,

Tom

1 Comment

Abed Alnaif

Abed Alnaif (view profile)

on 20 Mar 2015 at 14:51

Thanks, that's a neat solution!

tom_brot

tom_brot (view profile)

Contact us