Το SciPy (προφέρεται πίτα sai) είναι ένα πακέτο μαθηματικών βασισμένο σε numpy που περιλαμβάνει επίσης βιβλιοθήκες C και Fortran. Το SciPy μετατρέπει τη διαδραστική συνεδρία Python σε ένα πλήρες περιβάλλον επιστήμης δεδομένων όπως το MATLAB, το IDL, το Octave, το R ή το SciLab.
Σε αυτό το άρθρο, θα δούμε τις βασικές τεχνικές μαθηματικού προγραμματισμού - επίλυση προβλημάτων βελτιστοποίησης υπό όρους για μια βαθμωτή συνάρτηση πολλών μεταβλητών χρησιμοποιώντας το πακέτο scipy.optimize. Οι αλγόριθμοι βελτιστοποίησης χωρίς περιορισμούς έχουν ήδη συζητηθεί στο . Μπορείτε πάντα να λάβετε πιο λεπτομερή και ενημερωμένη βοήθεια για τις λειτουργίες scipy χρησιμοποιώντας την εντολή help(), Shift+Tab ή στο .
Εισαγωγή
Μια κοινή διεπαφή για την επίλυση προβλημάτων βελτιστοποίησης υπό όρους και χωρίς περιορισμούς στο πακέτο scipy.optimize παρέχεται από τη συνάρτηση minimize(). Ωστόσο, είναι γνωστό ότι δεν υπάρχει καθολική μέθοδος για την επίλυση όλων των προβλημάτων, επομένως η επιλογή μιας κατάλληλης μεθόδου, όπως πάντα, πέφτει στους ώμους του ερευνητή.
Ο κατάλληλος αλγόριθμος βελτιστοποίησης καθορίζεται χρησιμοποιώντας το όρισμα συνάρτησης minimize(..., method="").
Για τη βελτιστοποίηση υπό όρους μιας συνάρτησης πολλών μεταβλητών, είναι διαθέσιμες υλοποιήσεις των ακόλουθων μεθόδων:
trust-constr— αναζήτηση για ένα τοπικό ελάχιστο στην περιοχή εμπιστοσύνης. , ;SLSQP— διαδοχικός τετραγωνικός προγραμματισμός με περιορισμούς, Νευτώνεια μέθοδος για την επίλυση του συστήματος Lagrange. .TNC- Περικομμένος Newton Περιορισμένος, περιορισμένος αριθμός επαναλήψεων, κατάλληλος για μη γραμμικές συναρτήσεις με μεγάλο αριθμό ανεξάρτητων μεταβλητών. .L-BFGS-B— μια μέθοδος από την ομάδα Broyden–Fletcher–Goldfarb–Shanno, που υλοποιήθηκε με μειωμένη κατανάλωση μνήμης λόγω μερικής φόρτωσης διανυσμάτων από τον πίνακα Hessian. , .COBYLA— Περιορισμένη βελτιστοποίηση MARE με γραμμική προσέγγιση, περιορισμένη βελτιστοποίηση με γραμμική προσέγγιση (χωρίς υπολογισμό κλίσης). .
Ανάλογα με την επιλεγμένη μέθοδο, οι όροι και οι περιορισμοί για την επίλυση του προβλήματος τίθενται διαφορετικά:
- αντικείμενο τάξης
Boundsγια μεθόδους L-BFGS-B, TNC, SLSQP, trust-constr; - η λίστα
(min, max)για τις ίδιες μεθόδους L-BFGS-B, TNC, SLSQP, trust-constr; - ένα αντικείμενο ή μια λίστα αντικειμένων
LinearConstraint,NonlinearConstraintγια μεθόδους COBYLA, SLSQP, trust-constr. - λεξικό ή κατάλογος λεξικών
{'type':str, 'fun':callable, 'jac':callable,opt, 'args':sequence,opt}για μεθόδους COBYLA, SLSQP.
Περίληψη άρθρου:
1) Εξετάστε τη χρήση ενός αλγορίθμου βελτιστοποίησης υπό όρους στην περιοχή αξιοπιστίας (μέθοδος=”trust-constr”) με περιορισμούς που καθορίζονται ως αντικείμενα Bounds, LinearConstraint, NonlinearConstraint ;
2) Εξετάστε τον διαδοχικό προγραμματισμό χρησιμοποιώντας τη μέθοδο των ελαχίστων τετραγώνων (μέθοδος = "SLSQP") με περιορισμούς που καθορίζονται με τη μορφή λεξικού {'type', 'fun', 'jac', 'args'};
3) Αναλύστε ένα παράδειγμα βελτιστοποίησης κατασκευασμένων προϊόντων χρησιμοποιώντας το παράδειγμα ενός web studio.
Μέθοδος βελτιστοποίησης υπό όρους = "trust-constr"
Εφαρμογή της μεθόδου trust-constr βασισμένο στο για προβλήματα με περιορισμούς της μορφής της ισότητας και επί για προβλήματα με περιορισμούς με τη μορφή ανισοτήτων. Και οι δύο μέθοδοι υλοποιούνται από αλγόριθμους για την εύρεση ενός τοπικού ελάχιστου στην περιοχή εμπιστοσύνης και είναι κατάλληλες για προβλήματα μεγάλης κλίμακας.
Μαθηματική διατύπωση του προβλήματος της εύρεσης ενός ελάχιστου σε γενική μορφή:
Για αυστηρούς περιορισμούς ισότητας, το κάτω όριο ορίζεται ίσο με το άνω όριο .
Για περιορισμό μονής κατεύθυνσης, ορίζεται το ανώτερο ή το κατώτερο όριο np.inf με το αντίστοιχο πρόσημο.
Ας είναι απαραίτητο να βρούμε το ελάχιστο μιας γνωστής συνάρτησης Rosenbrock δύο μεταβλητών:
Σε αυτήν την περίπτωση, τίθενται οι ακόλουθοι περιορισμοί στον τομέα ορισμού του:
Στην περίπτωσή μας, υπάρχει μια μοναδική λύση στο σημείο , για τα οποία ισχύουν μόνο ο πρώτος και ο τέταρτος περιορισμός.
Ας δούμε τους περιορισμούς από κάτω προς τα πάνω και ας δούμε πώς μπορούμε να τους γράψουμε αναλυτικά.
Περιορισμοί и ας το ορίσουμε χρησιμοποιώντας το αντικείμενο Bounds.
from scipy.optimize import Bounds
bounds = Bounds ([0, -0.5], [1.0, 2.0])Περιορισμοί и Ας το γράψουμε σε γραμμική μορφή:
Ας ορίσουμε αυτούς τους περιορισμούς ως αντικείμενο LinearConstraint:
import numpy as np
from scipy.optimize import LinearConstraint
linear_constraint = LinearConstraint ([[1, 2], [2, 1]], [-np.inf, 1], [1, 1])Και τέλος ο μη γραμμικός περιορισμός σε μορφή πίνακα:
Ορίζουμε τον Jacobian matrix για αυτόν τον περιορισμό και έναν γραμμικό συνδυασμό του Hessian matrix με ένα αυθαίρετο διάνυσμα :
Τώρα μπορούμε να ορίσουμε έναν μη γραμμικό περιορισμό ως αντικείμενο NonlinearConstraint:
from scipy.optimize import NonlinearConstraint
def cons_f(x):
return [x[0]**2 + x[1], x[0]**2 - x[1]]
def cons_J(x):
return [[2*x[0], 1], [2*x[0], -1]]
def cons_H(x, v):
return v[0]*np.array([[2, 0], [0, 0]]) + v[1]*np.array([[2, 0], [0, 0]])
nonlinear_constraint = NonlinearConstraint(cons_f, -np.inf, 1, jac=cons_J, hess=cons_H)Εάν το μέγεθος είναι μεγάλο, οι πίνακες μπορούν επίσης να καθοριστούν σε αραιή μορφή:
from scipy.sparse import csc_matrix
def cons_H_sparse(x, v):
return v[0]*csc_matrix([[2, 0], [0, 0]]) + v[1]*csc_matrix([[2, 0], [0, 0]])
nonlinear_constraint = NonlinearConstraint(cons_f, -np.inf, 1,
jac=cons_J, hess=cons_H_sparse)ή ως αντικείμενο LinearOperator:
from scipy.sparse.linalg import LinearOperator
def cons_H_linear_operator(x, v):
def matvec(p):
return np.array([p[0]*2*(v[0]+v[1]), 0])
return LinearOperator((2, 2), matvec=matvec)
nonlinear_constraint = NonlinearConstraint(cons_f, -np.inf, 1,
jac=cons_J, hess=cons_H_linear_operator)Κατά τον υπολογισμό του Hessian matrix απαιτεί πολλή προσπάθεια, μπορείτε να χρησιμοποιήσετε μια τάξη . Οι ακόλουθες στρατηγικές είναι διαθέσιμες: BFGS и SR1.
from scipy.optimize import BFGS
nonlinear_constraint = NonlinearConstraint(cons_f, -np.inf, 1, jac=cons_J, hess=BFGS())Το Hessian μπορεί επίσης να υπολογιστεί χρησιμοποιώντας πεπερασμένες διαφορές:
nonlinear_constraint = NonlinearConstraint (cons_f, -np.inf, 1, jac = cons_J, hess = '2-point')Ο Jacobian πίνακας για περιορισμούς μπορεί επίσης να υπολογιστεί χρησιμοποιώντας πεπερασμένες διαφορές. Ωστόσο, σε αυτή την περίπτωση ο πίνακας της Έσσης δεν μπορεί να υπολογιστεί χρησιμοποιώντας πεπερασμένες διαφορές. Το Hessian πρέπει να οριστεί ως συνάρτηση ή χρησιμοποιώντας την κλάση HessianUpdateStrategy.
nonlinear_constraint = NonlinearConstraint (cons_f, -np.inf, 1, jac = '2-point', hess = BFGS ())Η λύση στο πρόβλημα βελτιστοποίησης μοιάζει με αυτό:
from scipy.optimize import minimize
from scipy.optimize import rosen, rosen_der, rosen_hess, rosen_hess_prod
x0 = np.array([0.5, 0])
res = minimize(rosen, x0, method='trust-constr', jac=rosen_der, hess=rosen_hess,
constraints=[linear_constraint, nonlinear_constraint],
options={'verbose': 1}, bounds=bounds)
print(res.x)`gtol` termination condition is satisfied.
Number of iterations: 12, function evaluations: 8, CG iterations: 7, optimality: 2.99e-09, constraint violation: 1.11e-16, execution time: 0.033 s.
[0.41494531 0.17010937]Εάν είναι απαραίτητο, η συνάρτηση για τον υπολογισμό του Hessian μπορεί να οριστεί χρησιμοποιώντας την κλάση LinearOperator
def rosen_hess_linop(x):
def matvec(p):
return rosen_hess_prod(x, p)
return LinearOperator((2, 2), matvec=matvec)
res = minimize(rosen, x0, method='trust-constr', jac=rosen_der, hess=rosen_hess_linop,
constraints=[linear_constraint, nonlinear_constraint],
options={'verbose': 1}, bounds=bounds)
print(res.x)ή το γινόμενο του Hessian και ενός αυθαίρετου διανύσματος μέσω της παραμέτρου hessp:
res = minimize(rosen, x0, method='trust-constr', jac=rosen_der, hessp=rosen_hess_prod,
constraints=[linear_constraint, nonlinear_constraint],
options={'verbose': 1}, bounds=bounds)
print(res.x)Εναλλακτικά, η πρώτη και η δεύτερη παράγωγος της βελτιστοποιούμενης συνάρτησης μπορούν να προσεγγιστούν. Για παράδειγμα, το Hessian μπορεί να προσεγγιστεί χρησιμοποιώντας τη συνάρτηση SR1 (οιονεί Νευτώνεια προσέγγιση). Η κλίση μπορεί να προσεγγιστεί με πεπερασμένες διαφορές.
from scipy.optimize import SR1
res = minimize(rosen, x0, method='trust-constr', jac="2-point", hess=SR1(),
constraints=[linear_constraint, nonlinear_constraint],
options={'verbose': 1}, bounds=bounds)
print(res.x)Μέθοδος βελτιστοποίησης υπό όρους = "SLSQP"
Η μέθοδος SLSQP έχει σχεδιαστεί για να επιλύει προβλήματα ελαχιστοποίησης μιας συνάρτησης με τη μορφή:
Πού и — σύνολα δεικτών εκφράσεων που περιγράφουν περιορισμούς με τη μορφή ισοτήτων ή ανισοτήτων. — σύνολα κάτω και άνω ορίων για τον τομέα ορισμού της συνάρτησης.
Οι γραμμικοί και μη γραμμικοί περιορισμοί περιγράφονται με τη μορφή λεξικών με κλειδιά type, fun и jac.
ineq_cons = {'type': 'ineq',
'fun': lambda x: np.array ([1 - x [0] - 2 * x [1],
1 - x [0] ** 2 - x [1],
1 - x [0] ** 2 + x [1]]),
'jac': lambda x: np.array ([[- 1.0, -2.0],
[-2 * x [0], -1.0],
[-2 * x [0], 1.0]])
}
eq_cons = {'type': 'eq',
'fun': lambda x: np.array ([2 * x [0] + x [1] - 1]),
'jac': lambda x: np.array ([2.0, 1.0])
}Η αναζήτηση για το ελάχιστο πραγματοποιείται ως εξής:
x0 = np.array([0.5, 0])
res = minimize(rosen, x0, method='SLSQP', jac=rosen_der,
constraints=[eq_cons, ineq_cons], options={'ftol': 1e-9, 'disp': True},
bounds=bounds)
print(res.x)Optimization terminated successfully. (Exit mode 0)
Current function value: 0.34271757499419825
Iterations: 4
Function evaluations: 5
Gradient evaluations: 4
[0.41494475 0.1701105 ]Παράδειγμα βελτιστοποίησης
Σε σχέση με τη μετάβαση στην πέμπτη τεχνολογική δομή, ας δούμε τη βελτιστοποίηση παραγωγής χρησιμοποιώντας το παράδειγμα ενός web studio, που μας φέρνει ένα μικρό αλλά σταθερό εισόδημα. Ας φανταστούμε τον εαυτό μας ως διευθυντή μιας γαλέρας που παράγει τρία είδη προϊόντων:
- x0 - πώληση σελίδων προορισμού, από 10 tr.
- x1 - εταιρικοί ιστότοποι, από 20 τρ.
- x2 - ηλεκτρονικά καταστήματα, από 30 τρ.
Η φιλική ομάδα εργασίας μας περιλαμβάνει τέσσερις junior, δύο μεσαίους και έναν ανώτερο. Το μηνιαίο ταμείο χρόνου εργασίας τους:
- Ιούνιος:
4 * 150 = 600 чел * час, - μεσαίες:
2 * 150 = 300 чел * час, - κύριος:
150 чел * час.
Αφήστε τον πρώτο διαθέσιμο junior να αφιερώσει (0, 1, 2) ώρες για την ανάπτυξη και την ανάπτυξη ενός ιστότοπου τύπου (x10, x20, x30), μεσαίου - (7, 15, 20), ανώτερου - (5, 10, 15 ) ώρες από την καλύτερη στιγμή της ζωής σας.
Όπως κάθε κανονικός διευθυντής, θέλουμε να μεγιστοποιήσουμε τα μηνιαία κέρδη. Το πρώτο βήμα προς την επιτυχία είναι να γράψετε την αντικειμενική συνάρτηση value ως το ποσό του εισοδήματος από προϊόντα που παράγονται ανά μήνα:
def value(x):
return - 10*x[0] - 20*x[1] - 30*x[2]Αυτό δεν είναι σφάλμα· κατά την αναζήτηση του μέγιστου, η αντικειμενική συνάρτηση ελαχιστοποιείται με το αντίθετο πρόσημο.
Το επόμενο βήμα είναι να απαγορεύσουμε στους υπαλλήλους μας να εργάζονται υπερβολικά και να εισάγουμε περιορισμούς στο ωράριο εργασίας:
Τι είναι ισοδύναμο:
ineq_cons = {'type': 'ineq',
'fun': lambda x: np.array ([600 - 10 * x [0] - 20 * x [1] - 30 * x[2],
300 - 7 * x [0] - 15 * x [1] - 20 * x[2],
150 - 5 * x [0] - 10 * x [1] - 15 * x[2]])
}Ένας επίσημος περιορισμός είναι ότι η παραγωγή προϊόντος πρέπει να είναι μόνο θετική:
bnds = Bounds ([0, 0, 0], [np.inf, np.inf, np.inf])Και, τέλος, η πιο ρόδινη υπόθεση είναι ότι λόγω της χαμηλής τιμής και της υψηλής ποιότητας, μια ουρά ικανοποιημένων πελατών δημιουργεί συνεχώς ουρά για εμάς. Μπορούμε να επιλέξουμε μόνοι μας τους μηνιαίους όγκους παραγωγής, με βάση την επίλυση του προβλήματος περιορισμένης βελτιστοποίησης scipy.optimize:
x0 = np.array([10, 10, 10])
res = minimize(value, x0, method='SLSQP', constraints=ineq_cons, bounds=bnds)
print(res.x)[7.85714286 5.71428571 3.57142857]Ας στρογγυλοποιήσουμε χαλαρά σε ακέραιους αριθμούς και ας υπολογίσουμε το μηνιαίο φορτίο των κωπηλατών με βέλτιστη κατανομή προϊόντων x = (8, 6, 3) :
- Ιούνιος:
8 * 10 + 6 * 20 + 3 * 30 = 290 чел * час; - μεσαίες:
8 * 7 + 6 * 15 + 3 * 20 = 206 чел * час; - κύριος:
8 * 5 + 6 * 10 + 3 * 15 = 145 чел * час.
Συμπέρασμα: για να λάβει ο σκηνοθέτης το άξιο του μέγιστο, είναι βέλτιστο να δημιουργείτε 8 σελίδες προορισμού, 6 ιστότοπους μεσαίου μεγέθους και 3 καταστήματα ανά μήνα. Σε αυτή την περίπτωση, ο ανώτερος πρέπει να οργώσει χωρίς να κοιτάξει ψηλά από το μηχάνημα, το φορτίο των μεσαίων θα είναι περίπου 2/3, οι junior λιγότερο από το μισό.
Συμπέρασμα
Το άρθρο περιγράφει τις βασικές τεχνικές για την εργασία με το πακέτο scipy.optimize, χρησιμοποιείται για την επίλυση προβλημάτων ελαχιστοποίησης υπό όρους. Προσωπικά χρησιμοποιώ scipy καθαρά για ακαδημαϊκούς σκοπούς, γι' αυτό και το παράδειγμα που δίνεται είναι τόσο κωμικού χαρακτήρα.
Πολλά θεωρητικά και εικονικά παραδείγματα μπορούν να βρεθούν, για παράδειγμα, στο βιβλίο του I.L. Akulich «Μαθηματικός προγραμματισμός σε παραδείγματα και προβλήματα». Περισσότερη σκληροπυρηνική εφαρμογή scipy.optimize για την κατασκευή μιας τρισδιάστατης δομής από ένα σύνολο εικόνων () μπορεί να προβληθεί σε .
Η κύρια πηγή πληροφοριών είναι όσοι επιθυμούν να συνεισφέρουν στη μετάφραση αυτής και άλλων ενοτήτων scipy Καλωσήρθες στο .
σας ευχαριστώ για συμμετοχή στην προετοιμασία της έκδοσης.
Πηγή: www.habr.com