Defining boundary condition for pde for pdepe function

I have got a pde as show in the function [c,f,s]. Unable to get the solution for the equation and not able figure out the mistake in boundary condition

L = 200;                 
s1 = 0.5; %equal to k at x=0               
s2 = 0;                  
T = 4;                   
qr = 0.218;               
f = 0.52;                
a = 0.0115;              
n = 2.03;                
ks = 31.6;               

x = linspace(0,L,100);
t = linspace(0,T,25);

options=odeset('RelTol',1e-4,'AbsTol',1e-4,'NormControl','off','InitialStep',1e-7)
u = pdepe(0,@unsatpde,@unsatic,@unsatbc,x,t,options,s1,s2,qr,f,a,n,ks);

I'M STRUCK AT THIS POINT. Following gives the editor .m files

% -------------------------------------------------------------------------
function [c,f,s] = trial1(x,t,u,DuDx)
global k n qr a ks p
c=1;
f = k*((DuDx)+1);
s = 0;
m =0.51;
q=qr+(p-qr)*(1+(-a*u)).^-m;
k=ks*((q-qr)/(p-qr))^0.5*(1-(1-((q-qr)/(p-qr))^(1/m))^m)^2;
% -------------------------------------------------------------------------
function u0 = unsatic(x,s1,s2,qr,f,a,n,ks)
u0 = 200+x; 

% -------------------------------------------------------------------------
function [pl,ql,pr,qr] = unsatbc(xl,ul,xr,ur,t,s1,s2,qr,f,a,n,ks)
pl = 0; 
ql = 1; 
pr = ur(1);
qr =0;

Answer 1

[Torsten]

1. "trial1" is not part of the list of functions you call pdepe with.

2. You will have to include s1,s2,qr,f,a,n,ks in the parameter list for "trial1".

3. in "trial1", you use k before you calculate it.

4. In unsatbc, you set as boundary conditions

u=0 at x=L

and

du/dx = -1/k at x=0

I don't know if this is what you want to set.

Best wishes

Torsten.

[vibha s]

Hi thanks for the suggestions,trying to extract the value of pL,qL and qR,pR in boundary condition code file. warning says error at intermediate points

[Torsten]

q will become complex because 1+(-a*u) becomes negative and is raised to the power of -m.

Best wishes

Torsten.

Answer 2

[Zana Taher]

Hi Torsten,

How would you make:

du/dx = 0 at x=0

instead of

du/dx = -1/k at x=0

?

[Torsten]

pl = -1.0;
ql = 1.0;

Answer 3

[jose luis huayanay villar]

uld you help me embed in the simulink? to carry out a control?