solving bvp differential equation

I want to solve bvp differential equation. but after a lot of effort I couldn't . my equation is:

(k1*cos(teta)^2+k3*sin(teta)^2)*(d^2(teta)/dz^2)+0.5*(k3-k1)*2*cos(teta)*sin(teta)*(d(teta)/dz)^2+(mu)^-1*kapa*b^2*sin(teta)*cos(teta)=0
NOTE:
1. k1, k3, kapa, mu are constants.
2. teta is a function of z.
3. my boundary condition is : teta(z=0)=0 ,  teta(z=d)=0 , teta'(z=d/2)=0

• I don't know to write third bondary condition ( my interval is [0,d]) *I try to solve it with the example of 3 in BVP_tutorial.pdf but I can't ( faced with a lot error)
• My code is:*

function mat4bvp
w = 2;
L=5*10^-6;   % is d
solinit = bvpinit(linspace(0,L,10),@mat4init,w);
sol = bvp4c(@mat4ode,@mat4bc,solinit);

fprintf('The fourth eigenvalue is approximately %7.3f.\n',...
      sol.parameters)

xint = linspace(0,L);
Sxint = deval(sol,xint);
plot(xint,Sxint(1,:))
axis([0 L 0 1.6])
title('Eigenfunction of Mathieu''s equation.') 
xlabel('x')
ylabel('solution y')
% -----------------------------------------------------------------------------
% w is magnetic field (B)
function dydx = mat4ode(x,y,w)
mu = 4*pi*10^-7;
kapa=9.5*10^-7;
k11=5.3*10^-12;
k33=7.3*10^-12;
dydx = [  y(2)
       -0.5*(k33-k11)*2*cos(y(x))*sin(y(x))*y(2)-((mu)^-1*kapa*w^2*sin(y(x))*cos(y(x)))];
% ------------------------------------------------------------------------------
function res = mat4bc(ya,yb,w)
res = [  ya(1) 
       yb(1)];
% ------------------------------------------------------------
function yinit = mat4init(x)
L=5*10^-6;
yinit = [cos((pi/L)*x)];