MATLAB: How to use pdepe with arbitrary initial condition

MATLABpdepe

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?

Best Answer

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);
Any code not stated here is identical to that in the question. It is trivial to turn this into a solve for input x,t vectors
*Matlab 2014b, I guess it will stop working at some point if pdepe.m gets upgraded.
Related Question