How to use pdepe with arbitrary initial condition

Question

Jonathan am 12 Aug. 2016

0
Verknüpfen

Direkter Link zu dieser Frage

https://de.mathworks.com/matlabcentral/answers/299432-how-to-use-pdepe-with-arbitrary-initial-condition

Kommentiert: Hafish Mahdi am 21 Nov. 2017

I want to use pdepe with an arbitrary initial condition, rather than a defined symbolic equation.

The problem is that pdepe does not allow this, due to some stuff around lines 229-239 in pdepe.m where it uses single value from the x-vector to calculate a single value of the initial condition rather than just taking the assigned initial condition wholesale.

The code I want would look something like this

function sol = actual_pde(init_c)
%  Solving u_t  = d((u_x)^2)/dx 
%               = 2*u_xx * u_x ;
%
% FYI: init_c = sin(x)
x = linspace(0,1,20);
t = linspace(0,.1,5);
sol = pdepe(0,@pd_pde,@pd_initc,@pd_bc,x,t);
% --------------------------------------------------------------
function [c,f,s] = pd_pde(x,t,u,DuDx)
c = 1;
f = abs(DuDx)*DuDx;
s = 0;
% --------------------------------------------------------------
function u0 = pd_initc(init_c)
 u0 = init_c
% --------------------------------------------------------------
function [pl,ql,pr,qr] = pd_bc(xl,ul,xr,ur,t)
pl = ul;
ql = 0;
pr = pi * exp(-t);
qr = 1;

but matlab will only accept

function u0 = pd_initc(x)
u0=sin(x);

It seems like the only solutions would be to edit the pdepe.m file around these lines:

% Form the initial values with a column of unknowns at each
% mesh point and then reshape into a column vector for dae.
temp = feval(pd_initc,xmesh(1),varargin{:});
if size(temp,1) ~= length(temp)
  error(message('MATLAB:pdepe:InvalidOutputICFUN'))
end
npde = length(temp);
y0 = zeros(npde,nx);
y0(:,1) = temp;
for j = 2:nx
  y0(:,j) = feval(pd_initc,xmesh(j),varargin{:});
end

But that seems like a huge workaround. Am I missing something obvious here?

0 Kommentare
-2 ältere Kommentare anzeigen-2 ältere Kommentare ausblenden

Melden Sie sich an, um zu kommentieren.

Melden Sie sich an, um diese Frage zu beantworten.

Answer 1

Jonathan am 15 Aug. 2016

0
Verknüpfen

Direkter Link zu dieser Antwort

https://de.mathworks.com/matlabcentral/answers/299432-how-to-use-pdepe-with-arbitrary-initial-condition#answer_231782

Bearbeitet: Jonathan am 15 Aug. 2016

In MATLAB Online öffnen

Well, here is a solution. It involves using global variables to get around matlab's implementation of the pdepe algorithm.

Obviously this is not best practice but, hey, it works*

function sol = pde_solve()
%  Solving u_t  = d((u_x)^2)/dx 
%               = 2*u_xx * u_x ;
%Setup vars
x = linspace(0,1,20);
t = linspace(0,.1,5);
%define global vars (Ugh) to get around limitations in defining initial
%conditions
global x_vector;
global init_c_vector
init_c_vector = sin(x);
x_vector = x;
% Do solve
sol = pdepe(0,@pd_pde,@pd_initc,@pd_bc,x,t);
% Extract the first solution component as u.
u = sol(:,:,1);
%%%%%%%%%%%%%%%%%%%
% Snipped out boring code
%%%%%%%%%%%%%%%%%%%
function u0 = pd_initc(x)
  %initialise global variables    
  global init_c_vector;
  global x_vector;
  % find the index of the x value pdepe wants to call
  [R,C] = find(x_vector==x,1,'first');
  %use the index to return the value of interest
  u0 = init_c_vector(R,C);