@samiyare You are right. That was my silly mistake!

@samiyare You can use continuation in N[bt]. You get results when N[bt] is 0.4. Then use continuation from 0.4 to 0.2, i.e. set
N[bt]:=cc*0.2+(1-cc)*0.4; #with cc in 0..1.
then N[bt] starts at cc=0 by being 0.4, where dsolve doesn't have a problem and works its way up to cc=1 where N[bt]=0.2. It happens to work.
You need only insert the optional argument 'continuation'=cc in dsolve.
I upload my version of your latest example.

MaplePrimes13-12-27o.mw

(Note the link should work now)

The error comes from dsolve. In your equation eq1 there is a term diff(u(eta),eta)/eta. You are considering the interval 0..1. Thus there will be a problem at eta = 0.
In dsolve you can specify the method as method=bvp[midrich]  (or middefer).
I didn't have any luck right away though. I got an error "Inintial Newton iteration not converging".
You probably need to do some exploration before embarking on the ambitious project of determining p2 and phi0 with the required properties. So just trying to make dsolve come out with some result for suitable values of p2 and phi0 would be a first step.
Comparing your equations with the equations considered on November 27:
http://www.mapleprimes.com/questions/200329-Automatic-Calculation-In-Dsolve-

I notice that there was no similar singularity. In that case the singularity would have been at eta=1 had it not been for the fact that the boundary condition was given at eta=1-zet, where then zet=0.5

Here you assign to zet twice. First zet:=0, secondly zet:=0.5: That of course just means that zet is the last value, but it does make you wonder. In Q the integrals are over the interval 0..1-zet. So that would be 0..0.5.

@Majmaj The output from print is NULL. That it prints to the screen does not in the usual Maple sense mean that it produces an output. You would have to copy and paste the result from the screen. You cannot save it to a variable.
You may want to refer to ?print where it states:
"The function print displays the values of the expressions appearing as arguments, and returns NULL as the function value. "

Try this:
restart;
generate_y := proc (y)
    y := (rand(1 .. 5))();
    print(y);
end proc;
qwerty:=12345;
generate_y(k);
%;
k;
yy:=generate_y(q);
#Compare with
generate_y2 := proc (y)
    y := (rand(1 .. 5))();
    NULL 
end proc;
zz:=generate_y2(p);

Just a follow-up.
The following extremely simple example shows an important difference between Maple 15 and Maple 16 and 17:

Maple 16 and 17:
B:=RootOf((x-n)*(x-2*n)*(x-3*n)*(x-4*n),x);
asympt(B,n,2);
                      
allvalues(%);
                       

Maple 15:
B:=RootOf((x-n)*(x-2*n)*(x-3*n)*(x-4*n),x);
asympt(B,n,2);
                        n

@John Fredsted At the danger of belaboring the point: The matrix M is only hermitian if the names in the entries are real.
op(M); #Showing shape=[]
type(M,'Matrix'(hermitian)); #Reasonable
LinearAlgebra:-HermitianTranspose(M) - M; #OK too
simplify(%) assuming real;
Only after assuming (as evalc does) that the entries are real do we get the zero matrix.

@Adri van der Meer You may want to accompany that by

M:=simplify(M) assuming real;

It would help a lot in getting people interested in your problem if they didn't have to enter the numbers themselves. You only provided an image.

@J4James If 8 is really a replacement for infinity, then one may speculate that the non-monotone solution of the ivp converges toward the monote solution on any bounded interval eta=0..E as "8" -> infinity.
If that is the case then the higher of the two f''(0) values would be the one to choose.

@Carl Love If so yes there are 2 solutions.
Try with 8 replaced by 12:
plot(p(.1,.1,.1,.1,f2),f2=-1..1,-.1..0.1); #Two zeros
r1:=fsolve(p(.1,.1,.1,.1,f2)=0,f2=-.3);
r2:=fsolve(p(.1,.1,.1,.1,f2)=0,f2=-.4);

@Carl Love Initially I just wanted to make use of the parameters option to see if that would save time. For that I had to turn the bvp into an ivp. Thus I needed f(0), f'(0), f''(0) given either concretely or as parameters.
The condition at eta=8 then had to be imposed afterwards.

@adel-00 I don't think we quite understand each other. Plotting beta versus t? Which t? I got the impression that that t was to be determined by the requirement (in this new case) N'(t) =0. 
What I was saying was that if for some t = t0 we had N'(t0) =0, then the right hand side of the ode (dsys1) would be zero, so either N(t0)=0 or N(t0) = a solution of a cubic (as in the previous problem).
So say we have N(t0) = N1 where N(t) = N1 makes the right hand side of the ode 0.
Then the ivp consisting of dsys1 and N(0)=N1 would have the constant solution N(t)=N1.
However, it also has the solution with initial condition N(0)=0.61299467707789. By uniqueness the solutions are the same, so we have that N1 =0.61299467707789. This determines beta. In this case we find:
subs(N(t)=rhs(ini1),rhs(dsys1));
solve(%=0,beta);
#returns beta roughly equal to zero.
My point, though, is that the problem as I understand it makes no sense.

@adel-00 You have a first order autonomous differential equation satisfying the conditions for uniqueness of the initial value problem. This implies that if y'(t) ever becomes 0 then it will have been zero at all times, i.e. that y is in fact constant.
For y' to be zero the right hand side of dsys1 has to be zero. The solutions to rhs = 0 are
y = 0 and the roots of the third order polynomial
pol:=-(1/3)*y^3+(1-epsilon)*y^2-(53/60)*y+(1-epsilon)/4
For some values of epsilon this polynomial has 3 real roots for other values only one.
You have y(0)=1/2 so if you want your solution to be constant that constant is 1/2.
Inserting y=1/2 into pol=0, you find that epsilon = 1/30.

@adel-00 This is an animation with epsilon as animation parameter. Each frame is a plot of y versus t.

This happens rather rarely and, as Markiyan points out, chances of getting an answer are greatly diminished. So I don't see any point in disallowing other languages.
However, I must admit to some irritation when somebody posts in a language other than English, in particular when it is in a language that I don't understand.