The syntax of eq1 is not correct You could use the construction `if`(condition, action1, action2) but I advice you to write insted

eq1 := piecewise(delta*(alpha+beta)*d*h >= x*d*h, Q2, x*d*h-delta*(alpha+beta)*d*h);

# more simply

eq1 := piecewise(delta*(alpha+beta) >= x, Q2, x-delta*(alpha+beta));

Whatever, eq1 is not an inequality (contrary to eq2..eq3a) and cannot be put between braces in Minimize command.

Other point eq3 is always verified as soon as eq3a is and is thus useless.

TR has only two indeterminates, I1 and I2, but inequalities eq2, and eq3a contain new indeterminates alpha, beta, d, h.
So it makes no sense include them into the Minimize command, just do

eq2a := I1 - I2 <= 0;
S := Minimize(TRC, {eq2a}, assume=nonnegative)

# or  Minimize(TRC, {eq2a}, I1=0..8, I2=0..8);  # assume=nonnegative is then useless
        S := S := [5.852170015000*10^12, [I1 = 0., I2 = 4.44089209850063*10^(-16)]]

Last point: I don't understand what you mean by "find what is Q1 and Q2...", please clarify.

Suppose Maple nor any other CAS does exist and just think a few seconds by yourself.

Some gives you this expression res which depends on 

indets(res, name)
            {k, lambda1, lambda2, n, omega, rho, t}

7 parameters and ask you to find the values of lambda1, lambda2, omega and rho such that numer(res)=0 and denom(res)<>0.

Note that they aretwwo extremely complex expressions:

length(numer(res)), length(denom(res))
                        67327531, 22899

How do you react?
Probably you would start saying  "What do you mean by solve with respect to lambda1, lambda2, omega and rho ? What about n, k, t ? Did you realize you have two equations and 4 unknowns?";
and next "Have you seen how complex those expressions are? How dare you ask for an exact solution?";
to end with "Your expression contains the imaginary unit, does it means some parameters are complex? What does denom(res)=0 mean: for instance has the real part to be 0 or has its module to be 0?"

Start defining exactly what you want... but do not ask for the impossible, as I said  Maple is not a crystal ball. 

I insist on the qualifier INITIAL as your questions seems to evolve according to the replies you receive (which is pretty annoying).

The attached file shows there is no problem at all (at least with Maple 2015.2) to fsolve the 6 equations and get positive real values:
fsolve_help_mmcdara.mw

Were you aware that the solutions verify lambda_d1 = lambda_d2 and  lambda_i1 = lambda_i2 ?

@MaPal93 


Trick.mw

@acer : it seems that DocumentTools is well adapted to place that "faked legend" you refered to

@Prakash J 

Shoot_Blasius-mmcdara_2.mw

(This is just an illustrative example for I don't know the values of the other parameters? So adjust my code to your own needs.)

plots[intersectplot] ?

Maybe using it will simplify your code

As X does not contain k it's hard to imagine how k disappears as you  "express X in terms of L".

@abdulganiy 

@abdulganiy 

I feel you are wrong: the value of h is not increasing, it remains equal to its initial value 0.01.
What is increasing is the product h*c which you fprint at column 1 under title "h".

You can verify this by writing

hhist    := Vector(round(N), h):

before the loop (hhist means "h history") and by rewritting the if structure where h is aimed to be changed:

 if Norm(err) <= tol then 
            h := 0.1 * h * (c + 1) * (tol/Norm(err))^(0.2);
            hhist[c] := h:
        end if;

Then verify after the loop the value of hhist...
In fact the condition Norm(err) <= tol is never met.

By the way, can you tell us what motivated you to write

h := 0.1 * h * (c + 1) * (tol/Norm(err))^(0.2)

Said this way it is a little bit obscure.
Can you provide something more explicit?

@janhardo 

Here is illustration where the value of printlevel is passed ato the procedure as an optional argument (default value = 1).
The file presents what happens in a "no-procedure" code and when this same code is encapsulated within a procedure.

In the "no-procedure" case you need printlevel:=3 to display all the intermediate operation but printlevel must be set to 8 in the "procedure" case. As I say at the end of the file I understand why we should set printlevel:=4, but I don't know where this value 8 comes from.

Nevertheless this illustrates how you could manage the display of intermediate results.
I don't know if this will help you.
Display_intermediate_results.mw

1. sol_f := dsolve({eq1, eq2, eq3, bcs, ...}, ....)
   is not sorrec: bcs is a set and dsolve requires aither a list or a set of expressions/equations.
   So you have to write sol_f := dsolve({eq1, eq2, eq3, op(bcs), ...}, ....) instead.
   0ther solution: define bcs this way bcs := {f(0) = 0, phi(0) = 1, theta(0) = 1, (D(f))(0) = 1}
    
2. Once corrected the execution of dsolve({eq1, eq2, eq3, op(bcs), ...}, numeric, ..., method = rkf45) you will get the explicit error error
   Error, (in dsolve/numeric/type_check) boundary value problem cannot be solved using the rkf45 method
   Remove method=rkf45 to get rid of it.
    
3. Once done another error appears
   Error, (in dsolve/numeric) 'parameters' must be specified as a list of unique unassigned names
   The reason is obvious: the option parameters is accepted only for ivp problems, not bvp problems.
   Solution: instanciate the system to solve with Parameters.
    

So let us do a simple thing whch encompasses points 1 to 3 above:

guess := [1$4];

bcs   := f(0) = 0, phi(0) = 1, theta(0) = 1, (D(f))(0) = 1:
sys   := eval({eq1, eq2, eq3, bcs}, Parameters):
sys   := sys union {f(100) = guess[1], (D(f))(100) = guess[2], theta(100) = guess[3], phi(100) = guess[4]}:

print~(sys): # You have 8 boundary conditions

                    ... odes not displayed, only are bcs...
                            f(0) = 0
                           f(100) = 1
                           phi(0) = 1
                          phi(100) = 1
                          theta(0) = 1
                         theta(100) = 1
                          D(f)(0) = 1
                         D(f)(100) = 1

# so the error: 
sol_f := dsolve(sys, numeric);
Error, (in dsolve/numeric/bvp/convertsys) too many boundary conditions: expected 7, got 8

Check if after having removed one bc (here the last one) dsolve runs without error:

sol_f := dsolve(sys[1..-2], numeric)
Error, (in dsolve/numeric/bvp) initial Newton iteration is not converging

Here appears another quite common problem. Execute help(dsolve, numeric_bvp, advanced) to see what this error means and find possible cures.
A classical cure often consists in introducing a continuation parameter. 

I can't help you more but I'm pretty sure some specialist here will.
In the meantime I advice you to think to this too many boundary conditions: expected 7, got 8 problem and to fix it yourself.

Good luck

Errors found

• In o_old :=Vector(n) n is undefined (is it a tyipo for N?)
   
• In the loop that follows: in for i from 1 to n n is indefined
   
• time is protected, you can use local time, but advice you chose another name of variable
   
• t0 is never initialized
   
• You also wrote something which appears like 
  Once rewritten in 1D input it means 
   0*_new[i] := o_old[i]+k[1, i]+2*k[2, i]+2*k[3, i]+k[4, i] 
  This produces the error
  Error, invalid left hand side of assignment
  So correct it (maybe you wanted to write o_new[i] := ... ?)
   
• ... it's likely I missed other mistakes or typos

I suggest you to modify your file accordingly to what I said and to come back to us after that.

Two last advices:

• As you are likely a newbie, do not use 2D mode and document style, prefer the classical 1D mode and worksheet. You will probably think it is old school but 1D worksheet enables avoiding the kind of typos your document seems to contain.
   
• Write loops in different blocks, for instance do this

  Block 1 = loop 1:
  Block 2 = loop 2:

  instead of 

  Block 1 = loop 1;
            loop 2:

  It is easier to check if the loops perform correctly with the first scheme, and if you get an error to find where it comes from.

@KIRAN SAJJAN 

In my answers it is clearly written that I removed the condition at Theta=-Pi/4.

So "My solution" verifies H(-Pi/4) = 1 and not (unless particular cases) H(Pi/4) = 1.
If you want to get H(Pi/4) = 1 remove the condition H(-Pi/4) = -1

Whatever, here are two files, read carefully the red highlighted text in both of them.
Complete_mmcdara_H_left_bc.mw
Complete_mmcdara_H_right_bc.mw

@KIRAN SAJJAN

Finally, I took pity on you and decided to do your job.
But the correct boundary conditions on H are up to you, and don't come asking me to help you with that, I'm done.

Look here Complete_mmcdara.mw