Il modo più efficiente per risolvere problemi di ottimizzazione parzialmente interi

Nov 04 2020

Permettere $0 < u < v < w < 10$essere numeri naturali. Voglio massimizzare$r$ in tutte le scelte di $u, v, w$ soggetto alle disuguaglianze $(42 - 100 r) - (6 - 10 r)(w + v) + vw(1 - r), -(42 - 100 r) + (6 - 10 r)(w + u) + uw(1 - r), (42 - 100 r) - (6 - 10 r)(u + v) + vu(1 - r) \ge 0.$ Per esempio, $(u,v,w)=(4,5,6)$ dovrebbe dare $r \le 4/15$ (questo è il primo $(u,v,w)$ Ho indovinato, ma non deve essere il migliore), e ora abbiamo bisogno $83$più valori. Senza la condizione dei numeri naturali, ho ricevuto un comando funzionante ma l'espressione risultante era brutta. Ho aggiunto le condizioni, ma ora il comando non è riuscito:

{
 {Maximize[{r, (42 - 100 r) - (6 - 10 r)*(w + v) + v*w *(1 - r) >= 
      0 && -(42 - 100 r) + (6 - 10 r)*(w + u) + u*w *(1 - r) >= 
      0 && (42 - 100 r) - (6 - 10 r)*(u + v) + v*u *(1 - r) >= 0, 
    0 < u < v < w < 10  , u \[Element] PositiveIntegers , 
    v \[Element] PositiveIntegers , w \[Element] PositiveIntegers  }, 
   r]},
 {\[Placeholder]}
}

C'è un modo per fare le cose in modo più efficiente in modo che il comando non fallisca? Se è possibile, vorrei anche sapere come restituire i valori di$u, v, w$ che raggiungono il massimo.

Ho appena provato

f[u_, v_, w_] := 
 Maximize[{r, (42 - 100 r) - (6 - 10 r)*(w + v) + v*w *(1 - r) >= 
     0 && -(42 - 100 r) + (6 - 10 r)*(w + u) + u*w *(1 - r) >= 
     0 && (42 - 100 r) - (6 - 10 r)*(u + v) + v*u *(1 - r) >= 0}, r]

in modo che potessi chiamare

Maximize[{f (u, v, w), 0 < u < v < w < 10, 
  u \[Element] PositiveIntegers, v \[Element] PositiveIntegers, 
  w \[Element] PositiveIntegers}, {u, v, w}]

successivo, ma quando ho definito f, ho ottenuto l'output $ Fallito.

Aggiornamento: ottengo

Maximize[{{4/15, {r -> 4/15}}, True, True, True, True}, {4, 5, 6}]

dopo aver definito la funzione e chiamato l'ultimo frammento di codice. Certamente, Mathematica non può avermi letto nel pensiero e sapere che ho già fatto il caso $ (4,5,6) $ e ho ottenuto $ 4/15, $ quindi sono tentato di credere che il primo valore che ho provato abbia davvero dare il massimo. Ma allora perché Mathematica restituisce \ $ Failed per ogni input? La dichiarazione di funzione aveva \$Failed, and this output had \$Fallito di lato.

Aggiornamento finale: in retrospettiva, quanto segue sarebbe stato più facile.

f[u, v, w] = Inverse[{{1, 1, 1}, {u, v, w}, {u^2, v^2, w^2}}]

NMaximize[{r, 
  VectorGreaterEqual[{f[u, v, w].{1 - r, 6 - 10 r, 42 - 100 r}, {0, 0,
       0}}] && 0 < u < v < w < 10 && u \[Element] PositiveIntegers && 
   v \[Element] PositiveIntegers && 
   w \[Element] PositiveIntegers}, {r, u, v, w}, 
 Method -> "DifferentialEvolution"]

Dopo aver ricevuto un errore 1/0 e "NMaximize: i seguenti vincoli non sono validi: ...", sostituire manualmente il vettore in VectorGuesdayEqual, rimuovere i denominatori dall'output ed eseguire di nuovo.

NMaximize[{r, {(42 - 100 r) (-v + w) + (6 - 10 r) (v^2 - w^2) + (1 - 
         r) (-v^2 w + v w^2), (42 - 100 r) (u - w) + (6 - 
         10 r) (-u^2 + w^2) + (1 - r) (u^2 w - u w^2), (42 - 
         100 r) (-u + v) + (6 - 10 r) (u^2 - v^2) + (1 - r) (-u^2 v + 
         u v^2)} \[VectorGreaterEqual] {0, 0, 0} && 
   0 < u < v < w < 10 && u \[Element] Integers && u > 0 && 
   v \[Element] Integers && v > 0 && w \[Element] Integers && 
   w > 0}, {r, u, v, w}, Method -> "DifferentialEvolution"] 

Risposte

3 user64494 Nov 04 2020 at 14:31

L'uso di NMinimizeinvece di Minimize

NMaximize[{r, (42 - 100 r) - (6 - 10 r)*(w + v) + v*w*(1 - r) >= 
0 && -(42 - 100 r) + (6 - 10 r)*(w + u) + u*w*(1 - r) >= 
0 && (42 - 100 r) - (6 - 10 r)*(u + v) + v*u*(1 - r) >= 0 && 
0 < u < v < w < 10 && u \[Element] PositiveIntegers && 
v \[Element] PositiveIntegers && w \[Element] PositiveIntegers}, {r, u, v, w}, 
Method -> "DifferentialEvolution"]
(*{0.321429, {r -> 0.321429, u -> 1, v -> 2, w -> 3}}*)

mentre

N[4/15]
(*0.266667*)

Aggiunta.

NMaximize[{r, (42 - 100 r) - (6 - 10 r)*(w + v) + v*w*(1 - r) >= 
0 && -(42 - 100 r) + (6 - 10 r)*(w + u) + u*w*(1 - r) >= 
0 && (42 - 100 r) - (6 - 10 r)*(u + v) + v*u*(1 - r) >= 0 && 
0 < u < v < w < 10 && u \[Element] PositiveIntegers && v \[Element] PositiveIntegers && 
w\[Element]PositiveIntegers},{r, u, v, w},Method-> {"DifferentialEvolution","ScalingFactor"->1}]
(*{0.95122, {r -> 0.95122, u -> 7, v -> 8, w -> 9}}*)

Molte grazie da parte mia a @joka per il prezioso avviso.

2 Akku14 Nov 04 2020 at 22:50

Ottieni una semplice soluzione analitica con l'aiuto di Reduce.

red = Reduce[(42 - 100 r) - (6 - 10 r)*(w + v) + v*w*(1 - r) >= 
 0 && -(42 - 100 r) + (6 - 10 r)*(w + u) + u*w*(1 - r) >= 
 0 && (42 - 100 r) - (6 - 10 r)*(u + v) + v*u*(1 - r) >= 0 && 
0 < u < v < w < 10 && u \[Element] Integers && 
v \[Element] Integers && w \[Element] Integers, r];

Maximize[{r, red}, {r, u, v, w}]

(*   {39/41, {r -> 39/41, u -> 7, v -> 8, w -> 9}}   *)
1 joka Nov 04 2020 at 18:04

Da $w,u,v$sono numeri interi e il loro intervallo comune ( rng-vedi codice) non è proibitivo (da 1 a 9, inclusi) è possibile scorrere tutte le combinazioni e selezionare quella migliore.

Il codice offerto di seguito, itera sui valori ammissibili per le variabili intere e registra una configurazione di variabili che rispettano i vincoli e massimizzano la funzione obiettivo se il valore ottimale è strettamente maggiore di qualsiasi valore registrato in precedenza.

(* common range of u, v, w *)
rng = Range[9]

(* formulating the constraints *)
poly1[r_] := 42 - 100 r
poly2[r_] := 6 - 10 r
poly3[r_] := 1 - r

c1[r_,u_,v_,w_] := poly1[r] - poly2[r] (w + v) + v w poly3[r]
c2[r_,u_,v_,w_] := -poly1[r] + poly2[r] (w + u) + u w poly3[r]
c3[r_,u_,v_,w_] := poly1[r] - poly2[r] (v + u) + u v poly3[r]

Preparare il loop

(* ranges for iterators *)
{i0, ie} = Through[{Min, Max}[rng]];
{j0, je} = {k0, ke} = {i0, ie};

(* store results in bag *)
vars = {u, v, w};
args = {0, 0, 0};

(* first entry is not a avalid sol *)
bag = {{0, Thread[vars -> args]}};

(* naive counter *)
c = 0;

Il ciclo effettivo; Monitorviene utilizzato per fornire un aggiornamento sull'iterazione corrente ( c).

Monitor[

 Do[

  ClearAll[c1N, c2N, c3N, cNs, objWCs, res, fMax, xOptRl, args];

  (* create constraints for given integer values of u,v,w  *)
  c1N[r_] := c1[r, i, j, k] // Evaluate;
  c2N[r_] := c2[r, i, j, k] // Evaluate;
  c3N[r_] := c3[r, i, j, k] // Evaluate;

  (* prep the inequality constraints *)
  cNs[r_] := {c1N[r], c2N[r], c3N[r]};

  (* list of objective with all constraints *)
  objWCs[r_] := Join[{r}, Thread[cNs[r] >= 0], {i < j, j < k}] // Evaluate;

  (* the actual optimization under constraints *)
  res = Check[NMaximize[objWCs[r], r], $Failed] // Quiet; Which[ (* there was a result *) res // FailureQ /* Not, ( {fMax, xOptRl} = res; If[ (* new opt better than current opt *) fMax > (bag // Last /* First), args = Thread[vars -> {i, j, k}]; (* update res *) bag = Flatten[{bag, {{fMax, args}}}, 1] ] ), (* no opt *) True, $Failed

   ];

  c += 1, {i, i0, ie}, {j, j0, je}, {k, k0, ke}], c]

Dopo aver valutato il ciclo sopra, la riga seguente

bag // Rest 

valuta a

{
  {0.321429, {u -> 1, v -> 2, w -> 3}},
  {0.333333, {u -> 5, v -> 7, w -> 8}},
  {0.4, {u -> 5, v -> 8, w -> 9}},
  {0.5, {u -> 6, v -> 7, w -> 8}},
  {0.75, {u -> 6,v -> 8,w -> 9}},
  {0.95122, {u -> 7, v -> 8, w -> 9}}
 } 

Ovviamente, il valore ottimale per rsotto i vincoli è 0.95122e questo può essere ottenuto quando u -> 7, v -> 8e w -> 9.

Modifica per accogliere una considerazione minore sollevata nei commenti

Per affrontare il problema sollevato da @ user64494 nei commenti, è possibile, ma non necessario , sostituirlo NMaximizecon Maximizeinside the loop e ottenere il valore ottimale di rin forma razionale.

Effettuando la sostituzione e valutando, produce

 {
  {9/28, {u -> 1, v -> 2, w -> 3}}, 
  {1/3, {u -> 5, v -> 7, w -> 8}}, 
  {2/5, {u -> 5, v -> 8, w -> 9}}, 
  {1/2, {u -> 6, v -> 7, w -> 8}}, 
  {3/4, {u -> 6, v -> 8, w -> 9}}, 
  {39/41, {u -> 7, v -> 8, w -> 9}}
 }

Ora, valutando le seguenti righe di codice,

(* make a list of rules for the values of r, u, v, w at the optimum *)
sol = bag // Last /* (MapAt[Rule[r, #] &, #, 1] &) /* Flatten;

(* make a list of the non-negativity constraint functions *)
cs = Map[RightComposition[Apply[Sequence], #] &, {c1, c2, c3}];

(* evaluate the non-negativity constraints at the optimum solution *)
Through[cs[{r, u, v, w} /. sol]] >= 0 // Thread /* Apply[And]

produce True.

Ciò significa che tutti i vincoli di non negatività sono soddisfatti alla soluzione ottima ottenuta sopra.

Ovviamente, il u < v < wvincolo è anche banalmente adattato alla soluzione ottenuta.

Semplice controllo dei valori di $u,v,w$al livello ottimale rivela che in effetti 7 < 8 < 9è anche True .