Skip to content

Commit

Permalink
ruff
Browse files Browse the repository at this point in the history
  • Loading branch information
@Paul-Saves
Paul-Saves committed Dec 17, 2024
1 parent cfbd735 commit 4d2e8d7
Showing 1 changed file with 92 additions and 57 deletions.
149 changes: 92 additions & 57 deletions tutorial/SBO/SMT_EGO_application.ipynb
View file
Original file line number Diff line number Diff line change
Expand Up @@ -35808,11 +35808,13 @@
"goal = Branin(ndim=ndim)\n",
"constraint = Sphere(ndim=ndim)\n",
"\n",
"x0_lim = [-5,10]\n",
"x1_lim = [0,15]\n",
"x0_lim = [-5, 10]\n",
"x1_lim = [0, 15]\n",
"\n",
"xlimits = np.array([x0_lim, x1_lim])\n",
"sampling = LHS(xlimits=xlimits, random_state = 12) #Random state is the seed number, only used for reproducibility \n",
"sampling = LHS(\n",
" xlimits=xlimits, random_state=12\n",
") # Random state is the seed number, only used for reproducibility\n",
"num = 15\n",
"x_t = sampling(num)\n",
"\n",
Expand Down Expand Up @@ -35871,6 +35873,7 @@
],
"source": [
"from smt.surrogate_models import KRG\n",
"\n",
"# Training the surrogate for the objective function\n",
"y_hat = KRG(theta0=[1e-2])\n",
"y_hat.set_training_values(x_t, y_t)\n",
Expand Down Expand Up @@ -35901,37 +35904,46 @@
"source": [
"# Create plotting data\n",
"num = 500\n",
"x = np.zeros([num,2])\n",
"x = np.zeros([num, 2])\n",
"x[:, 0] = np.linspace(-5.0, 10.0, num)\n",
"x[:, 1] = np.linspace(0.0, 15.0, num)\n",
"\n",
"X, Y = np.meshgrid(x[:, 0], x[:, 1]) # Create a grid\n",
"f = np.zeros([num,num])\n",
"g = np.zeros([num,num])\n",
"X, Y = np.meshgrid(x[:, 0], x[:, 1]) # Create a grid\n",
"f = np.zeros([num, num])\n",
"g = np.zeros([num, num])\n",
"for i in range(num):\n",
" f[:,i] = y_hat.predict_values(np.array([X[:,i], Y[:,i]]).T).flatten()\n",
" g[:,i] = g_hat.predict_values(np.array([X[:,i], Y[:,i]]).T).flatten()\n",
" f[:, i] = y_hat.predict_values(np.array([X[:, i], Y[:, i]]).T).flatten()\n",
" g[:, i] = g_hat.predict_values(np.array([X[:, i], Y[:, i]]).T).flatten()\n",
"\n",
"from IPython.display import clear_output\n",
"\n",
"clear_output()\n",
"fig, axes = plt.subplots(1, 2, figsize=(12, 4))\n",
"# Ploting objective function\n",
"contour1 = axes[0].contourf(X, Y, f, levels=40, cmap='RdBu_r') # Filled contours in red and blue\n",
"contour_lines = axes[0].contour(X, Y, f, levels=40, colors='black') # Contour lines in black\n",
"contour1 = axes[0].contourf(\n",
" X, Y, f, levels=40, cmap=\"RdBu_r\"\n",
") # Filled contours in red and blue\n",
"contour_lines = axes[0].contour(\n",
" X, Y, f, levels=40, colors=\"black\"\n",
") # Contour lines in black\n",
"fig.colorbar(contour1, ax=axes[0])\n",
"axes[0].scatter(x_t[:, 0], x_t[:, 1], color = 'black')\n",
"axes[0].set_title('Predicted object function', fontsize=16)\n",
"axes[0].set_xlabel('x1')\n",
"axes[0].set_ylabel('x2')\n",
"axes[0].scatter(x_t[:, 0], x_t[:, 1], color=\"black\")\n",
"axes[0].set_title(\"Predicted object function\", fontsize=16)\n",
"axes[0].set_xlabel(\"x1\")\n",
"axes[0].set_ylabel(\"x2\")\n",
"\n",
"# Ploting constraint function\n",
"contour2 = axes[1].contourf(X, Y, g, levels=40, cmap='RdBu_r') # Filled contours in red and blue\n",
"contour_lines = axes[1].contour(X, Y, g, levels=40, colors='black') # Contour lines in black\n",
"contour2 = axes[1].contourf(\n",
" X, Y, g, levels=40, cmap=\"RdBu_r\"\n",
") # Filled contours in red and blue\n",
"contour_lines = axes[1].contour(\n",
" X, Y, g, levels=40, colors=\"black\"\n",
") # Contour lines in black\n",
"fig.colorbar(contour2, ax=axes[1])\n",
"axes[1].scatter(x_t[:, 0], x_t[:, 1], color = 'black')\n",
"axes[1].set_title('Predicted constraint function', fontsize=16)\n",
"axes[1].set_xlabel('x1')\n",
"axes[1].set_ylabel('x2')"
"axes[1].scatter(x_t[:, 0], x_t[:, 1], color=\"black\")\n",
"axes[1].set_title(\"Predicted constraint function\", fontsize=16)\n",
"axes[1].set_xlabel(\"x1\")\n",
"axes[1].set_ylabel(\"x2\")"
]
},
{
Expand Down Expand Up @@ -35962,18 +35974,24 @@
"lambda_i = 40\n",
"y_min = np.min(y_t[np.where(g_t <= lambda_i)[0]])\n",
"\n",
"\n",
"def ExpectedImprovement(x, objective, y_min):\n",
" from scipy.special import erf\n",
"\n",
" mu = objective.predict_values(x)\n",
" var = objective.predict_variances(x)+1e-6\n",
" EI = (y_min - mu)*(0.5+0.5*erf((y_min - mu) / np.sqrt(var * 2))) + np.sqrt(var/(2*np.pi))*np.exp(-pow(y_min-mu,2)/(2*var))\n",
" var = objective.predict_variances(x) + 1e-6\n",
" EI = (y_min - mu) * (0.5 + 0.5 * erf((y_min - mu) / np.sqrt(var * 2))) + np.sqrt(\n",
" var / (2 * np.pi)\n",
" ) * np.exp(-pow(y_min - mu, 2) / (2 * var))\n",
" return EI\n",
"\n",
"\n",
"def ProbabilityOfFeasibility(x, constraint, lambda_i):\n",
" from scipy.special import erf\n",
"\n",
" mu = constraint.predict_values(x)\n",
" var = constraint.predict_variances(x)+1e-6\n",
" PF = (0.5+0.5*erf((lambda_i - mu) / np.sqrt(var * 2)))\n",
" var = constraint.predict_variances(x) + 1e-6\n",
" PF = 0.5 + 0.5 * erf((lambda_i - mu) / np.sqrt(var * 2))\n",
" return PF"
]
},
Expand Down Expand Up @@ -36004,14 +36022,15 @@
],
"source": [
"from pymoo.algorithms.soo.nonconvex.es import ES\n",
"from pymoo.problems import get_problem\n",
"from pymoo.optimize import minimize\n",
"\n",
"cEI_problem = lambda x: ExpectedImprovement(x, y_hat, y_min)*ProbabilityOfFeasibility(x, g_hat, lambda_i)\n",
"def cEI_problem(x):\n",
" return ExpectedImprovement(x, y_hat, y_min) * ProbabilityOfFeasibility(x, g_hat, lambda_i)\n",
"\n",
"import numpy as np\n",
"from pymoo.core.problem import Problem\n",
"\n",
"\n",
"class MyProblem(Problem):\n",
" def __init__(self, object, constraint, y_min, lambda_i):\n",
" self.object = object\n",
Expand All @@ -36021,17 +36040,16 @@
" super().__init__(n_var=2, n_obj=1, xl=np.array([-5, 0]), xu=np.array([10, 15]))\n",
"\n",
" def _evaluate(self, x, out, *args, **kwargs):\n",
" cEI = ExpectedImprovement(x, self.object, self.y_min)*ProbabilityOfFeasibility(x, self.constraint, self.lambda_i)\n",
" cEI = ExpectedImprovement(\n",
" x, self.object, self.y_min\n",
" ) * ProbabilityOfFeasibility(x, self.constraint, self.lambda_i)\n",
" out[\"F\"] = -cEI\n",
"\n",
"\n",
"problem = MyProblem(y_hat, g_hat, y_min, lambda_i)\n",
"algorithm = ES(n_offsprings=200, rule=1.0 / 7.0)\n",
"\n",
"res = minimize(problem,\n",
" algorithm,\n",
" (\"n_gen\", 200),\n",
" seed=1,\n",
" verbose=False)\n",
"res = minimize(problem, algorithm, (\"n_gen\", 200), seed=1, verbose=False)\n",
"\n",
"clear_output()\n",
"print(\"Best solution found: \\nX = %s\\nF = %s\" % (res.X, res.F))"
Expand Down Expand Up @@ -36062,45 +36080,62 @@
],
"source": [
"# Create plotting data\n",
"EI_M = np.zeros([num,num])\n",
"PF_M = np.zeros([num,num])\n",
"cEI_M = np.zeros([num,num])\n",
"EI_M = np.zeros([num, num])\n",
"PF_M = np.zeros([num, num])\n",
"cEI_M = np.zeros([num, num])\n",
"for i in range(num):\n",
" EI = ExpectedImprovement(np.array([X[:,i], Y[:,i]]).T, y_hat, y_min)\n",
" PF = ProbabilityOfFeasibility(np.array([X[:,i], Y[:,i]]).T, g_hat, lambda_i)\n",
" EI_M[:,i] = EI.flatten()+1e-12 # The 1e-12 is only added to ensure a pretty plot and does not affect the result\n",
" PF_M[:,i] = PF.flatten()\n",
" cEI_M[:,i] =(EI*PF).flatten()+1e-12 # The 1e-12 is only added to ensure a pretty plot and does not affect the result\n",
" EI = ExpectedImprovement(np.array([X[:, i], Y[:, i]]).T, y_hat, y_min)\n",
" PF = ProbabilityOfFeasibility(np.array([X[:, i], Y[:, i]]).T, g_hat, lambda_i)\n",
" EI_M[:, i] = (\n",
" EI.flatten() + 1e-12\n",
" ) # The 1e-12 is only added to ensure a pretty plot and does not affect the result\n",
" PF_M[:, i] = PF.flatten()\n",
" cEI_M[:, i] = (\n",
" (EI * PF).flatten() + 1e-12\n",
" ) # The 1e-12 is only added to ensure a pretty plot and does not affect the result\n",
"\n",
"from IPython.display import clear_output\n",
"\n",
"n_levels = 20\n",
"clear_output()\n",
"fig, axes = plt.subplots(1, 2, figsize=(12, 4))\n",
"# Ploting objective function\n",
"contour1 = axes[0].contourf(X, Y, EI_M, levels=n_levels, cmap='RdBu_r') # Filled contours in red and blue\n",
"contour_lines = axes[0].contour(X, Y, EI_M, levels=n_levels, colors='black') # Contour lines in black\n",
"contour1 = axes[0].contourf(\n",
" X, Y, EI_M, levels=n_levels, cmap=\"RdBu_r\"\n",
") # Filled contours in red and blue\n",
"contour_lines = axes[0].contour(\n",
" X, Y, EI_M, levels=n_levels, colors=\"black\"\n",
") # Contour lines in black\n",
"fig.colorbar(contour1, ax=axes[0])\n",
"axes[0].set_title('Expected Improvement', fontsize=16)\n",
"axes[0].set_xlabel('x1')\n",
"axes[0].set_ylabel('x2')\n",
"axes[0].set_title(\"Expected Improvement\", fontsize=16)\n",
"axes[0].set_xlabel(\"x1\")\n",
"axes[0].set_ylabel(\"x2\")\n",
"\n",
"# Ploting constraint function\n",
"contour2 = axes[1].contourf(X, Y, PF_M, levels=n_levels, cmap='RdBu_r') # Filled contours in red and blue\n",
"contour_lines = axes[1].contour(X, Y, PF_M, levels=n_levels, colors='black') # Contour lines in black\n",
"contour2 = axes[1].contourf(\n",
" X, Y, PF_M, levels=n_levels, cmap=\"RdBu_r\"\n",
") # Filled contours in red and blue\n",
"contour_lines = axes[1].contour(\n",
" X, Y, PF_M, levels=n_levels, colors=\"black\"\n",
") # Contour lines in black\n",
"fig.colorbar(contour2, ax=axes[1])\n",
"axes[1].set_title('Probability of Feasibility', fontsize=16)\n",
"axes[1].set_xlabel('x1')\n",
"axes[1].set_ylabel('x2')\n",
"axes[1].set_title(\"Probability of Feasibility\", fontsize=16)\n",
"axes[1].set_xlabel(\"x1\")\n",
"axes[1].set_ylabel(\"x2\")\n",
"\n",
"fig, axes = plt.subplots(figsize=(6, 4))\n",
"# Ploting objective function\n",
"contour1 = axes.contourf(X, Y, cEI_M, levels=n_levels, cmap='RdBu_r') # Filled contours in red and blue\n",
"contour_lines = axes.contour(X, Y, cEI_M, levels=n_levels, colors='black') # Contour lines in black\n",
"axes.scatter([res.X[0]], [res.X[1]], marker='*', color = 'yellow', s = 50)\n",
"contour1 = axes.contourf(\n",
" X, Y, cEI_M, levels=n_levels, cmap=\"RdBu_r\"\n",
") # Filled contours in red and blue\n",
"contour_lines = axes.contour(\n",
" X, Y, cEI_M, levels=n_levels, colors=\"black\"\n",
") # Contour lines in black\n",
"axes.scatter([res.X[0]], [res.X[1]], marker=\"*\", color=\"yellow\", s=50)\n",
"fig.colorbar(contour1, ax=axes)\n",
"axes.set_title('constrained Expected Improvement', fontsize=16)\n",
"axes.set_xlabel('x1')\n",
"axes.set_ylabel('x2')"
"axes.set_title(\"constrained Expected Improvement\", fontsize=16)\n",
"axes.set_xlabel(\"x1\")\n",
"axes.set_ylabel(\"x2\")"
]
}
],
Expand Down

0 comments on commit 4d2e8d7

Please sign in to comment.