Question: Contour plot pf ode (couple equations)

Is it possible the contour plot of these coupled equations

restart:with(plots): Nb := .2; Nt := .2; h1 := 1+a*sin(x);h2:= -b-d*sin(x); a := 1;b:=0.1;d:=0.1;F:=1-Q-d;Q:=-1;beta:=0.1;alpha:=0.1;G:=0.5;B:=0.5; eq1 := diff(alpha*f(y),y,y,y,y)+G*(diff(theta(y),y,y))+ B*(diff(phi(y),y,y))+6* beta*(diff(f(y),y,y))*(diff(f(y),y,y,y))^(2)+3* beta*(diff(f(y),y,y,y,y))*(diff(f(y),y,y))^(2)=0; eq2 := diff(theta(y), y, y)+Nb*(diff(theta(y), y))*(diff(phi(y), y))+Nt*(diff(theta(y), y))^2 = 0; eq3 := diff(phi(y), y, y)+Nb*(diff(theta(y), y, y))/Nt = 0;

Boundary conditions

bc := f(h1) = (1/2)*F, f(h2) = -(1/2)*F, (D(f))(h1) = -1, theta(h2) = 1, phi(h2) = 1, (D(f))(h2) = -1, theta(h1) = 0, phi(h1) = 0; f(1 + sin(x)) = 0.9500000000, f(-0.1 - 0.1 sin(x)) = -0.9500000000, D(f)(1 + sin(x)) = -1, theta(-0.1 - 0.1 sin(x)) = 1, phi(-0.1 - 0.1 sin(x)) = 1, D(f)(-0.1 - 0.1 sin(x)) = -1, theta(1 + sin(x)) = 0, phi(1 + sin(x)) = 0 M[x,y] := dsolve({eq1,eq2,eq3,bc}, numeric);

p1 := contourplot(M[x,y], [x,y, f(y)],x=1..10,.4;y=1..10,.4)); plots[display](p1);