@Andiguys 

This can be easily proved by examining the Hessian of Pi.
It has a null eigenvalue and the condition for tit to be (semi)negative definiteis never realized with your data.
So there is no maximum inside the domain D = (0, 1)x(0, 1)x(0, 1) of variation of (Pn, Pr, w) and Pi reaches its maximum value on some point of the boundary of D.

More of this one can prove that a critical point in D is always a saddle point whatever the data you use.
So, to be clear: you will never get a macimum of Pi within D whatever your data (unless I did something wrong in my demonstration)

Proposal_follow_up_sand15.mw

@JoyDivisionMan

I see you didn't even take time to say if my answer satisfied you.
Is that a simple omission, or rudeness?

@Andiguys 

And whart are the variables relarive to which you want to minimize the objective function?
What about the others: can you provide numeric values or them?

@Andiguys 

Could you simply post a worksheet wihich contains only the objective function and the constraints instead of mixing your new problem with what I did previously?

By the way: it seems quite common for you to ask a question and, once you have received an answer, to ask a completely different question without even saying whether the solution we gave you suits you or not.
I find this particularly unpleasant.

@vv 

And to be clear I am not pleading for transparency or nominative voting.
Besides, it is not unusual for one member to say to another, “I vote up” and that seems sufficient to me.

@mthkvv 

How or where did you find this expression (MMA maybe)?

@Ronan @vv

Your solutions ar far better than mine, I vote up for both of you.

@Andiguys 

Point 1 "I am simultaneously solving three equations with three variables"
There is no reason to do this unless you want to make more complicated the job of Maple.
You already solved K2 versus w: why don't you plug this solution into the system {K3, K4} and solve a sysem of two equations in two unknowns?

Point 2 "This process is giving me three different solution sets"
No, in fact you have 4 solutions and not 3: one of them contains RootOfs which, after having applied allvalues, provides 2 equations.
So you have 4 solutions in total.

Point 3 "I am unsure which one is correct or should be considered"
Who do you think could be sure of the correct solution to consider if you are not capable to be yourself?

Point 4 "All parameters are positive and greater than zero. All decision variable Pn, Pr and w should be positive equal to or greater than zero"
Consider the extremely simple case of the solution of K2 wrt w. It writes

w = (Pr+i2-tau-Crm ) / 2

Assuming Pr > 0, i2 > 0, tau > 0, Crm > 0 does not imply that w >= 0.
So your positivity constraints on "All the parameters" is completely useless for choosing the "correct" (as you say) solution.
More of this you implicitely assume there is only one "correct" solution Pn >= 0, Pr >= 0, w >= 0, but the attached file proves you can get 0, 1, 2 "correct" solutions and maybe even 3 or 4.

new_Q_solve_sand15.mw

@salim-barzani 

Read again the file I sent you in this COMMENT

If it is not clear enough use this one Line-1-Done_sand15.mw : this is exactly the same file but the corrections are highlighted this way: 𝛌, r  instead of  𝛌₁, 𝛌₂, r₁ i!!!

@salim-barzani 

Login:
I used to use my personnal account (mmcdara) when I was home and my professional one (sand15) at the office.
But since the major Mapleprimes outage this summer  mmcdara is no longer accessible and I fell back on sand15 only.

Error
The function U is not correctly defined, see the simple preliminary example and the ensuing correct definition of U in this worksheet: Line-1-Done_sand15.mw

Some equations contain things like (D[1](f))(x, y, 0, t) because you never told Maple what f(x, y, z, t) was. So it cannot evalluate this D term which annoys you.

You should avoid writting things like 

eqf := f(x, y, z, t) = theta1*theta2+Bij(1, 2):

eq17 := u(x, y, z, t) = 2*(diff(f(x, y, z, t), x))/f(x, y, z, t):

z := 0:
ans := solve({diff(rhs(equ), x), diff(rhs(equ), y)}, {x, y}, explicit)

and replace them by 

f := unapply(theta1*theta2+Bij(1, 2), (x, y, z, t)):             

u := unapply(2*(diff(f(x, y, z, t), x))/f(x, y, z, t), (x, y, z, t)): 

du := diff~(u(x, y, 0, t), [x, y]):   # which simply compute derivatives at z=0 (z remains z)
ans := solve(du, [x, y], explicit)

While doing this you truly define functions f and u and so things like  (D[1](f))(x, y, 0, t) are automatically evaluated and bever appear.

Here is the beginning of what I would do (just the beginning because at some point I no longer understand what you want to achieve)

eval_sand15.mw

By the way, you might be relieved to know that the Error, (in plots:-display) expecting plot structure but received: contour2 doesn't exist with maple 2015: eval_Maple_2015.mw

@DEE_Engineering 

I'm glad I could be of help to you.

Here is the quite general procedure FitRationalFraction which fits any rational fraction P(x)/Q(x) where P(x) and Q(x) are dense polynomial functions with respective degrees p and q.
This "has-to-be-dense" limitation comes from the choice of defining this rational fraction in a very concise way (indeed the list [p, q] is passed to FitRationalFraction) and could be easily removed.

Four examples are provided, the last one uses your data whose result is



FitRationalFraction.mw

@Preben Alsholm 

Should this be considered as an example of code regression or do you consider it's a marginal drawback of recent versions?

sol := dsolve(simplify(ode, trig)):

But doing this returns two solutions (additionally of lengthy expressions)...

ode.mw

@Preben Alsholm 

I thought to make a reply to @alice abut was unable to reproduce the issue with Maple 2015:

restart
kernelopts(version);
           Maple 2015.2, APPLE UNIVERSAL OSX, Dec 20 2015, Build ID 1097895
Digits;
                               10

int(1e6*exp(-1e6*t), t=0..infinity);
                               1.

And for even larger lambda values

int(1e60*exp(-1e60*t), t=0..infinity);
                               1.