Advanced Search

Andiguys

100 Reputation

5 Badges

1 years, 260 days
Contact Andiguys

MaplePrimes Activity


These are questions asked by Andiguys

Multiple curves in a dual-axis plot ...

Andiguys 100 Maple 2019
dualaxisplot plotting
November 09 2025
2 8

I would like to combine all the plots into a single figure. The curves S1, S2, and S3 represent the manufacturer’s profit as Ce​ varies, and S12, S22, and S33 represent the retailer’s profit for the same changes in Ce​. I want all of these displayed together in one plot using a dual y-axis: one axis for the manufacturer’s profit and the other for the retailer’s profit, with Ce on the x-axis. How to create such a dual-axis plot with appropriate scaling so that the differences between the curves are also clearly visible.

> 

restart
> 

with(Optimization); with(plots); with(Student[VectorCalculus]); with(LinearAlgebra)
> 

``
> 

_local(Pi)

Pi

 (1)
> 

`Π_12` := (0.1455251030e-2*Ce+.5352049476)*(0.369876310e-1-0.3638127575e-2*Ce)+(.8*(-.1671790360+1.121361872*Ce))*(0.1849381518e-1-0.1819063782e-2*Ce)-Ce*(0.1849381518e-1-0.1819063782e-2*Ce)

(0.1455251030e-2*Ce+.5352049476)*(0.369876310e-1-0.3638127575e-2*Ce)+(-.1337432288+.8970894976*Ce)*(0.1849381518e-1-0.1819063782e-2*Ce)-Ce*(0.1849381518e-1-0.1819063782e-2*Ce)

 (2)
> 

`Π_22` := (0.1455251030e-2*Ce+.5356096675)*(0.355258312e-1-0.3638127575e-2*Ce)+(.8*(-.1184158360+1.121361872*Ce))*(0.1776291535e-1-0.1819063782e-2*Ce)-Ce*(0.1776291535e-1-0.1819063782e-2*Ce)

(0.1455251030e-2*Ce+.5356096675)*(0.355258312e-1-0.3638127575e-2*Ce)+(-0.9473266880e-1+.8970894976*Ce)*(0.1776291535e-1-0.1819063782e-2*Ce)-Ce*(0.1776291535e-1-0.1819063782e-2*Ce)

 (3)
> 

`Π_32` := (0.1455251030e-2*Ce+.5356038465)*(0.355403838e-1-0.3638127575e-2*Ce)+(.8*(-.1179012835+1.121361872*Ce))*(0.1777019161e-1-0.1819063782e-2*Ce)-Ce*(0.1777019161e-1-0.1819063782e-2*Ce)

(0.1455251030e-2*Ce+.5356038465)*(0.355403838e-1-0.3638127575e-2*Ce)+(-0.9432102680e-1+.8970894976*Ce)*(0.1777019161e-1-0.1819063782e-2*Ce)-Ce*(0.1777019161e-1-0.1819063782e-2*Ce)

 (4)
> 

S12 := plot(`Π_12`, Ce = 0 .. 0.9e-1, color = [red], labels = ["Ce", "Manufacturer Profit"], labeldirections = ["horizontal", "vertical"], legend = [`#msubsup(mi("Pi"),mi("m"),mn("W"));`]); S22 := plot(`Π_22`, Ce = 0 .. 0.9e-1, color = [green], labels = ["Ce", "Manufacturer profit"], labeldirections = ["horizontal", "vertical"], legend = [`#msubsup(mi("Pi"),mi("m"),mn("D"));`]); S32 := plot(`Π_32`, Ce = 0 .. 0.9e-1, color = [blue], labels = ["Ce", "Manufacturer profit"], labeldirections = ["horizontal", "vertical"], legend = [`#msubsup(mi("Pi"),mi("m"),mn("S"));`])

 

 

 
> 

`Π_1` := (-0.60726413e-1*Ce+.6173851967)*(0.1849381518e-1-0.1819063782e-2*Ce)-0.2500000000e-1*(0.1849381518e-1-0.1819063782e-2*Ce)^2

(-0.60726413e-1*Ce+.6173851967)*(0.1849381518e-1-0.1819063782e-2*Ce)-0.2500000000e-1*(0.1849381518e-1-0.1819063782e-2*Ce)^2

 (5)
> 

`Π_2` := (-0.60726413e-1*Ce+.5929853242)*(0.1776291535e-1-0.1819063782e-2*Ce)-0.2500000000e-1*(0.1776291535e-1-0.1819063782e-2*Ce)^2

(-0.60726413e-1*Ce+.5929853242)*(0.1776291535e-1-0.1819063782e-2*Ce)-0.2500000000e-1*(0.1776291535e-1-0.1819063782e-2*Ce)^2

 (6)
> 

`Π_3` := (-0.60726413e-1*Ce+.5932282299)*(0.1777019161e-1-0.1819063782e-2*Ce)-0.2500000000e-1*(0.1777019161e-1-0.1819063782e-2*Ce)^2

(-0.60726413e-1*Ce+.5932282299)*(0.1777019161e-1-0.1819063782e-2*Ce)-0.2500000000e-1*(0.1777019161e-1-0.1819063782e-2*Ce)^2

 (7)
> 

S1 := plot(`Π_1`, Ce = 0 .. 0.9e-1, color = [yellow], labels = ["Ce", "Retailer profit"], labeldirections = ["horizontal", "vertical"], legend = [`#msubsup(mi("Pi"),mi("r"),mn("W"));`]); S2 := plot(`Π_2`, Ce = 0 .. 0.9e-1, color = [black], labels = ["Ce", "Retailer profit"], labeldirections = ["horizontal", "vertical"], legend = [`#msubsup(mi("Pi"),mi("r"),mn("D"));`]); S3 := plot(`Π_3`, Ce = 0 .. 0.9e-1, color = [grey], labels = ["Ce", "Retailer profit"], labeldirections = ["horizontal", "vertical"], legend = [`#msubsup(mi("Pi"),mi("r"),mn("S"));`])

 

 

 
> 

dualaxisplot(plot(`Π_22`, Ce = 0 .. 0.9e-1, color = ["red"], labels = ["Ce", "Manufacturer profit"], labeldirections = ["horizontal", "vertical"], legend = [`#msubsup(mi("Pi"),mi("m"),mn("D"));`]), plot(`Π_2`, Ce = 0 .. 0.9e-1, color = ["green"], labels = ["Ce", "Retailer profit"], labeldirections = ["horizontal", "vertical"], legend = [`#msubsup(mi("Pi"),mi("r"),mn("D"));`]), title = "fairnes cost Comparison")

 
> 

display({S1, S12, S2, S22, S3, S32})

 
 

``

Download All_plots_Combined.mw

Isolating positive and negative terms in...

Andiguys 100 Maple 2019
optimization solve inequality
November 07 2025
1 2

How can I isolate the positive terms on one side and the negative terms on the other? Is there a systematic way to split and rearrange the expression so that I can determine the conditions under which L is positive or negative?

> 

restart
> 

``
> 

L := simplify(-rho0*(-Cn*Cr*alpha*d*rho0+Cr*Pr*alpha*d*rho0-Cr*alpha*d*delta*rho0+Ce*alpha*d*delta+Cr*alpha*d*rho0-Crm*alpha*d*delta+Pr*alpha*d*delta-alpha*c*d*delta-alpha*delta*g*i2-Ce*alpha*d-Ce*d*delta-2*Cn*alpha*rho0+Crm*alpha*d-Pr*alpha*d+2*Pr*alpha*rho0+2*a*alpha*delta+alpha*c*d-alpha*c*delta-2*alpha*delta*rho0+alpha*g*i2+Ce*d-2*a*alpha+alpha*c+2*alpha*rho0)/(Cr*alpha*d*delta*rho0^2-Cr*alpha*d*rho0^2+2*alpha*delta*rho0^2-2*alpha*rho0^2-d*delta^2+2*d*delta-d)) = 0

-2*((((-(1/2)*rho0*Cr-(1/2)*Crm-(1/2)*c+(1/2)*Ce+(1/2)*Pr)*delta-(1/2)*Cr*(Cn-Pr-1)*rho0+(1/2)*Crm+(1/2)*c-(1/2)*Ce-(1/2)*Pr)*d+(-(1/2)*i2*g+a-rho0-(1/2)*c)*delta+(-Cn+Pr+1)*rho0+(1/2)*i2*g-a+(1/2)*c)*alpha-(1/2)*d*Ce*(delta-1))*rho0/((rho0^2*(Cr*d+2)*alpha-d*(delta-1))*(delta-1)) = 0

 (1)
 

``

Download Positive_Negative_Isolation.mw

if statement inside a do-loop...

Andiguys 100 Maple 2019
coordinates
October 31 2025
1 5

How can I use an if statement inside a do - loop in maple and then plot a figure based on the conditional results? I need help with the correct syntax for combining the do-loop, `if` condition, and plotting command.

Sheet attached below:
Algorithmm.mw

μ₁ and μ₂ conditions in constrained opti...

Andiguys 100 Maple 2019
optimization lagrange-multipliers duplicate-question
September 29 2025
1 10

I am optimizing an objective function subject to two constraints. I have obtained the solution, but I am unsure how to determine the conditions for the Lagrange multipliers μ1​ and μ2​. For example, using the KKT conditions, if μ1>0 then the solution is X, and if μ2>0 then the solution is Y. Could you guide me on how to find μ1 and μ2 along with their corresponding conditions?

> 

restart
> 

with(Optimization); with(plots); with(LinearAlgebra)
> 

_local(Pi)

Pi

 (1)
 
> 

`π_m` := proc (i1) options operator, arrow; (w-i1)*(1/2+(1/2)*(i1-i2)/tau)*(1-(Pn-Pr)/(1-delta))-C0-(1/2)*Cr*(1/2+(1/2)*(i1-i2)/tau)^2*(1-(Pn-Pr)/(1-delta))^2+Ce*rho0*(1-(Pn-Pr)/(1-delta)) end proc
> 

C1 := (2*rho0-1)*tau+i2 <= i1

(2*rho0-1)*tau+i2 <= i1

 (2)
> 

C2 := i1 <= tau+i2

i1 <= tau+i2

 (3)
> 

NULL
> 

NULL
> 

``
> 

NULL
> 

# No equality constraints
#
# Inequality constraints must be of the form f[i](i1) <= 0
#
# Thus
#
# (1) For C1
#     so
      f[1] := (i1) -> (2*rho0 - 1)*tau + i2-i1:
      
#
# (2) C2
#     so
      f[2] := (i1) -> i1-(tau + i2):
      
#
> 

# Lagrangian (we want to maximize `&pi;_m` so to minimize -`&pi;_m`

L := -`&pi;_m`(i1) + add(f[i](i1)*mu[i], i=1..2):
> 

dLdi1 := collect(diff(L, i1), [i1]):
> 

KKT_conditions := [
                    seq(mu[i] >= 0, i=1..2),         # Dual feasibility conditions
                    dLdi1 = 0,                       # Stationarity condition
                    seq(``(f[i](i1)) <= 0, i=1..2),  # Primal feasibility conditions
                    add(mu[i]*f[i](i1) = 0, i=1..2)  # Complementary slackness
                  ]:

print~(KKT_conditions);

0 <= mu[1]

  

0 <= mu[2]

  

((1-(Pn-Pr)/(1-delta))/tau+(1/4)*Cr*(1-(Pn-Pr)/(1-delta))^2/tau^2)*i1+(1/2-(1/2)*i2/tau)*(1-(Pn-Pr)/(1-delta))-(1/2)*w*(1-(Pn-Pr)/(1-delta))/tau+(1/2)*Cr*(1/2-(1/2)*i2/tau)*(1-(Pn-Pr)/(1-delta))^2/tau-mu[1]+mu[2] = 0

  

``((2*rho0-1)*tau+i2-i1) <= 0

  

``(i1-tau-i2) <= 0

  

((2*rho0-1)*tau+i2-i1)*mu[1]+(i1-tau-i2)*mu[2] = 0

  

[]

 (4)


Analysis of dLdp1

> 

with(LargeExpressions):

DLDi1 := collect(dLdi1, i1, Veil[K]);

(1/4)*K[1]*i1-(1/4)*K[2]

 (5)
> 

# Thus dLdi1 = 0 verifies

isolate((5), i1)
 

i1 = K[2]/K[1]

 (6)
> 

# Where K1 and K2 are

i := 'i':
KS := [ seq(K[i] = Unveil[K](K[i]), i=1..LastUsed[K] ) ];

[K[1] = (-1+delta+Pn-Pr)*(Cr*Pn-Cr*Pr+Cr*delta+4*delta*tau-Cr-4*tau)/(tau^2*(-1+delta)^2), K[2] = (4*delta^2*tau^2*mu[1]-4*delta^2*tau^2*mu[2]+Cr*Pn^2*i2-Cr*Pn^2*tau-2*Cr*Pn*Pr*i2+2*Cr*Pn*Pr*tau+2*Cr*Pn*delta*i2-2*Cr*Pn*delta*tau+Cr*Pr^2*i2-Cr*Pr^2*tau-2*Cr*Pr*delta*i2+2*Cr*Pr*delta*tau+Cr*delta^2*i2-Cr*delta^2*tau+2*Pn*delta*i2*tau-2*Pn*delta*tau^2+2*Pn*delta*tau*w-2*Pr*delta*i2*tau+2*Pr*delta*tau^2-2*Pr*delta*tau*w+2*delta^2*i2*tau-2*delta^2*tau^2+2*delta^2*tau*w-8*delta*tau^2*mu[1]+8*delta*tau^2*mu[2]-2*Cr*Pn*i2+2*Cr*Pn*tau+2*Cr*Pr*i2-2*Cr*Pr*tau-2*Cr*delta*i2+2*Cr*delta*tau-2*Pn*i2*tau+2*Pn*tau^2-2*Pn*tau*w+2*Pr*i2*tau-2*Pr*tau^2+2*Pr*tau*w-4*delta*i2*tau+4*delta*tau^2-4*delta*tau*w+4*tau^2*mu[1]-4*tau^2*mu[2]+Cr*i2-Cr*tau+2*i2*tau-2*tau^2+2*tau*w)/(tau^2*(-1+delta)^2)]

 (7)


Analysis of the two constraints

> 

beta

beta

 (8)
> 

Cs := { seq(beta[i] = subs(i1=0, f[i](i1)), i=1..2) };

{beta[1] = (2*rho0-1)*tau+i2, beta[2] = -tau-i2}

 (9)
> 

 

> 

# two constraints problem

Cmin_value := min(rhs~(Cs));  # the "true" constraints
Cmax_value := max(rhs~(Cs)):

Cmin_name  := min(lhs~(Cs));  # the "abstract" constraints
Cmax_name  := max(lhs~(Cs)):

Cmin_name  := LowerBound;  # and an even more abstract form
Cmax_name  := UpperBound:

Complementary_Slackness := mu[1]*(Cmin_name-i1) , mu[2]*(i1-Cmax_name);

L2     := -`&pi;_m`(i1) + add(Complementary_Slackness):
dL2di1 := diff(L2, i1):
dL2di1 := map(simplify, %, size);
 

min(-tau-i2, (2*rho0-1)*tau+i2)

  

min(beta[1], beta[2])

  

LowerBound

  

mu[1]*(LowerBound-i1), mu[2]*(i1-UpperBound)

  

(1/2)*(tau+i1-i2)*(-1+delta+Pn-Pr)/(tau*(-1+delta))+(1/2)*(-w+i1)*(-1+delta+Pn-Pr)/(tau*(-1+delta))+(1/4)*Cr*(tau+i1-i2)*(-1+delta+Pn-Pr)^2/(tau^2*(-1+delta)^2)-mu[1]+mu[2]

 (10)
> 

type(dL2di1, linear([i1, seq(mu[i], i=1..2)]));

SC_CS_sols := solve(
  {
    dL2di1 = 0,
    op([Complementary_Slackness] =~ 0)
  }
  ,
  {i1, seq(mu[i], i=1..2)}
):

true

  

Main: Entering solver with 3 equations in 3 variables
Main: attempting to solve as a linear system
Main: attempting to solve as a polynomial system
Main: Polynomial solver successful. Exiting solver returning 1 solution

  
> 

# How many solutions did we get?

numelems([SC_CS_sols]);

# And those solutions are charecterized by

map(s -> if eval(lambda[1], s) = 0 then
           if eval(lambda[2], s) = 0 then
             "i1 belongs to interval (LowerBound, UpperBound)"
           else
             "i1 is equal to the UpperBound"
           end if
         else
           "i1 is equal to the LowerBound"
         end if
         , [SC_CS_sols]);

# So we get as expected the three possible situations.

3

  

["i1 is equal to the LowerBound", "i1 is equal to the LowerBound", "i1 is equal to the LowerBound"]

 (11)
> 

Cmin_value

min(-tau-i2, (2*rho0-1)*tau+i2)

 (12)
 

 

Download I_1_Optimum_condition.mw

Taking so much time to solve?...

Andiguys 100 Maple 2019
solve syntax parametric
September 26 2025
2 17

Easy way to solve 3 equation simultaneously?:

> 

restart
> 

with(Optimization); with(plots); with(LinearAlgebra)
 
> 

K2 := -(Pn-Pr)/(1-delta)+(-beta*((-2*delta+2)*tau^2+((-Pn+Pr-delta+1)*Cr+(-2+2*delta)*(w+i2))*tau+Cr*i2*(-1+delta+Pn-Pr))*upsilon/((4*delta-4)*tau+Cr*(-1+delta+Pn-Pr))+Pr)/delta+(Pr-w-Cr)*beta*(-2+2*delta)*tau*upsilon/(((4*delta-4)*tau+Cr*(-1+delta+Pn-Pr))*delta) = 0

 

 (1)
> 

K3 := 1-(Pn-Pr)/(1-delta)-(Pn-Cn)/(1-delta)+(Pr-w-Cr)*(1/(1-delta)-(-beta*(Cr*i2-Cr*tau)*upsilon/((4*delta-4)*tau+Cr*(-1+delta+Pn-Pr))+beta*((-2*delta+2)*tau^2+((-Pn+Pr-delta+1)*Cr+(-2+2*delta)*(w+i2))*tau+Cr*i2*(-1+delta+Pn-Pr))*upsilon*Cr/((4*delta-4)*tau+Cr*(-1+delta+Pn-Pr))^2)/delta)+Ce*rho0/(1-delta) = 0

 

 (2)
> 

K4 := (Pn-Cn)/(1-delta)+(Pn-Pr)/(1-delta)-(-beta*((-2*delta+2)*tau^2+((-Pn+Pr-delta+1)*Cr+(-2+2*delta)*(w+i2))*tau+Cr*i2*(-1+delta+Pn-Pr))*upsilon/((4*delta-4)*tau+Cr*(-1+delta+Pn-Pr))+Pr)/delta+(Pr-w-Cr)*(-1/(1-delta)-(-beta*(-Cr*i2+Cr*tau)*upsilon/((4*delta-4)*tau+Cr*(-1+delta+Pn-Pr))-beta*((-2*delta+2)*tau^2+((-Pn+Pr-delta+1)*Cr+(-2+2*delta)*(w+i2))*tau+Cr*i2*(-1+delta+Pn-Pr))*upsilon*Cr/((4*delta-4)*tau+Cr*(-1+delta+Pn-Pr))^2+1)/delta)-Ce*rho0/(1-delta) = 0

 

 (3)
> 

sol := solve({K2, K3, K4}, {Pn, Pr, w})

Download Q_solve.mw

1 2 3 4 5 6 7 Last Page 3 of 11