Experiment No. 9: Simulated Annealing (SA) for Engineering Optimization
Course: PEML3001 – Decision Making and Optimization Laboratory
Target Level: M.Tech (Project Engineering & Management / Mechanical Engineering)
1. Aim & Objectives
Aim
To design, implement, and analyze the Simulated Annealing (SA) optimization algorithm for solving a constrained non-linear engineering optimization problem (Pressure Vessel Cost Minimization) and compare its performance with conventional gradient-based methods and population-based metaheuristics such as Genetic Algorithms (GA) and Particle Swarm Optimization (PSO).
Learning Objectives
After completing this experiment, you will be able to:
- Understand Physical Analogies: Map thermodynamic annealing concepts (thermodynamic state, temperature, free energy) directly to optimization parameters (candidate solution, control parameter, cost function).
- Handle Probabilistic Search: Apply the Metropolis criterion to balance exploration (global search) and exploitation (local refinement).
- Overcome Local Minima: Formulate controlled randomization strategies that allow hill-climbing out of non-convex local traps.
- Formulate Constrained Design Problems: Model engineering design variables (continuous, integer, and discrete) along with complex stress and geometric constraints using exterior dynamic penalty functions.
- Analyze Parameter Sensitivity: Assess how initial temperature T_0, cooling rate \alpha, and neighborhood step size dictate convergence rates and solution quality.
- Compare Metaheuristics: Critically evaluate performance tradeoffs between trajectory-based metaheuristics (SA) and population-based algorithms (GA, PSO).
2. Theory & Physical Metaphor
Simulated Annealing (SA) is a trajectory-based stochastic metaheuristic proposed by Kirkpatrick, Gelatt, and Vecchi (1983) and independently by Cerny (1985). It adapts the Metropolis-Hastings algorithm (Metropolis et al., 1953) from statistical mechanics to general optimization problems.
The Metallurgical Analogy
In physical annealing, a solid material is heated above its melting point to a high temperature where atoms move freely in a high-energy disordered state. The material is then cooled extremely slowly (annealed).
- Slow Cooling: Allows atoms sufficient time to redistribute and arrange themselves into a minimum energy, highly ordered crystalline lattice.
- Rapid Cooling (Quenching): Traps thermal fluctuations, leading to a metastable, high-energy amorphous or brittle state with internal stresses (analogous to getting trapped in a local optimum).
Physical Annealing (Metallurgy) Optimization Domain (Simulated Annealing)
---------------------------------- -----------------------------------------
Material State / Atom Configuration <---> Candidate Solution Vector (x)
Internal Energy (E) <---> Objective / Cost Function f(x)
Temperature (T) <---> Control Parameter / Temperature (T)
Low Energy Lattice State <---> Global Optimal Solution (x*)
Quenching (Rapid Cooling) <---> Greedy / Local Descent Trajectory
3. Mathematical Background
3.1 Mathematical Formulation of Constrained Optimization
A general non-linear constrained engineering optimization problem is stated as:
To handle constraints within SA, an Exterior Penalty Function transforms the problem into an unconstrained objective \Phi(\mathbf{x}, T):
where \lambda, \mu \gg 0 are large penalty multipliers.
3.2 Acceptance Probability (The Metropolis Criterion)
Let \mathbf{x}{\text{current}} be the current solution state with cost E{\text{current}} = f(\mathbf{x}{\text{current}}), and \mathbf{x}{\text{new}} be a perturbed candidate solution with cost E_{\text{new}} = f(\mathbf{x}{\text{new}}). The change in system energy is:
The probability P of accepting the candidate solution \mathbf{x}{\text{new}} is governed by:
To decide whether an uphill move (\Delta E > 0) is accepted, a uniformly distributed random number r \sim U(0, 1) is generated:
Acceptance Probability P(ΔE, T) vs. ΔE
1.0 |------------------
| \ High T (Exploration: Accepts large positive ΔE)
| \ Medium T
| \ Low T (Exploitation: Approaching local greedy search)
0.0 +--------------------------> ΔE (Worse solution cost increase)
3.3 Cooling Schedules
The cooling schedule specifies how the control parameter T drops over time k.
| Cooling Schedule | Mathematical Formula | Key Characteristics |
|---|---|---|
| Geometric (Exponential) | T_{k+1} = \alpha \cdot T_k \quad (\alpha \in [0.80, 0.99]) | Most popular; easy to implement; predictable runtime. |
| Linear | T_{k+1} = T_k - \eta | Rapid cooling at low T; prone to premature convergence. |
| Logarithmic (Lundy-Mees) | T_{k+1} = \dfrac{T_k}{1 + \beta T_k} | Asymptotically guarantees convergence to global optimum; very slow. |
| Adaptive (Feedback) | T_{k+1} = T_k \cdot \left(1 - \dfrac{\gamma \cdot \sigma_E}{T_k}\right) | Adjusts cooling rate dynamically based on standard deviation of energy \sigma_E. |
4. Benchmark Problem: Pressure Vessel Cost Minimization
The algorithm is applied to the ASME Pressure Vessel Design Benchmark (Kannada & Turner, 1991). The objective is to minimize total fabrication, material, and welding costs of a cylindrical vessel with hemispherical heads.
+------------------------------------+
/| |\
/ | | \
| | <-------------- L ---------------> | |
\ | | /
\| |/
+------------------------------------+
|<- Ts R ->| Th ->| |
Design Variables
- x_1 = T_s: Shell thickness (integer multiple of 0.0625\text{ in.})
- x_2 = T_h: Head thickness (integer multiple of 0.0625\text{ in.})
- x_3 = R: Inner radius (10.0 \le R \le 200.0\text{ in.}, continuous)
- x_4 = L: Length of cylindrical section (10.0 \le L \le 200.0\text{ in.}, continuous)
Objective Function
Constraints
-
-
-
-
5. Algorithm & Flowchart
-
-
-
Step-by-Step Execution Sequence
- Initialization: Set initial parameters: T_0, T_{\min}, cooling rate \alpha, epoch length N_{\text{epoch}}. Select an initial feasible solution vector \mathbf{x}^{(0)} and compute E_0 = \Phi(\mathbf{x}^{(0)}). Set \mathbf{x}_{\text{best}} = \mathbf{x}^{(0)}.
- Outer Loop: While T > T_{\min}:
- Inner Loop (Markov Chain): Perform N_{\text{epoch}} iterations:
- Perturbation: Generate neighbor \mathbf{x}' = \mathbf{x} + \boldsymbol{\delta}, where \boldsymbol{\delta} \sim U(-\boldsymbol{\Delta}, \boldsymbol{\Delta}).
- Boundary Enforcers: Enforce integer grid bounds for discrete variables (x_1, x_2) and clamp continuous bounds (x_3, x_4).
- Evaluation: Compute new energy E' = \Phi(\mathbf{x}') and delta \Delta E = E' - E_{\text{current}}.
- Metropolis Decision:
- If \Delta E \le 0 \implies \text{Accept } \mathbf{x}'.
- If \Delta E > 0 \implies \text{Generate } r \sim U(0,1). If r < \exp(-\Delta E / T), accept \mathbf{x}'; otherwise retain current state.
- Global Best Update: If \Phi(\mathbf{x}{\text{accepted}}) < \Phi(\mathbf{x}{\text{best}}), set \mathbf{x}{\text{best}} = \mathbf{x}{\text{accepted}}.
- Cooling Step: Update temperature: T = \alpha \cdot T.
- Termination: Output global optimal decision vector \mathbf{x}{\text{best}} and minimum cost f(\mathbf{x}{\text{best}}).
[ START ]
|
[ Initialize: x0, T = T0, ]
[ T_min, alpha, N_epoch ]
|
+-------------->|
| (Is T > T_min?)
| |
| +-----+-----+
| YES NO
| | |
| [ Iteration = 0 ] +---> [ Output x_best, f(x_best) ]
| | |
| +---->| [ STOP ]
| | (Iter < N_epoch?)
| | |
| | +--+--+
| | YES NO -------------------+
| | | |
| | [ Generate Neighbor x' ] |
| | [ Compute ΔE = f(x')-f(x)] |
| | | |
| | (Is ΔE <= 0?) |
| | | |
| |+--+--+ |
| |YES NO |
| || | |
| || (r < exp(-ΔE/T)?) |
| || | |
| || +--+--+ |
| || YES NO |
| || | | |
| |v v v |
| [Accept] [Reject] |
| | |
| [Update x_best if better] |
| | |
| [Iter = Iter + 1] |
| | |
+---+ v
[ T = alpha * T ]
6. Implementation (Python 3.x)
import numpy as np
import math
class SimulatedAnnealingPressureVessel:
def __init__(self, T0=10000.0, Tmin=1e-3, alpha=0.95, N_epoch=100):
self.T0 = T0
self.Tmin = Tmin
self.alpha = alpha
self.N_epoch = N_epoch
# Variable Bounds: x1 (Ts), x2 (Th), x3 (R), x4 (L)
self.bounds = [(1, 99), (1, 99), (10.0, 200.0), (10.0, 200.0)]
def objective_cost(self, x):
# x1 and x2 are discrete multiples of 0.0625
x1 = x[0] * 0.0625
x2 = x[1] * 0.0625
x3 = x[2]
x4 = x[3]
cost = (0.6224 * x1 * x3 * x4 +
1.7781 * x2 * (x3**2) +
3.1661 * (x1**2) * x4 +
19.84 * (x1**2) * x3)
return cost
def constraints_penalty(self, x):
x1 = x[0] * 0.0625
x2 = x[1] * 0.0625
x3 = x[2]
x4 = x[3]
g1 = -x1 + 0.0193 * x3
g2 = -x2 + 0.00954 * x3
g3 = -np.pi * (x3**2) * x4 - (4.0/3.0) * np.pi * (x3**3) + 1296000.0
g4 = x4 - 240.0
pen = 0.0
for g in [g1, g2, g3, g4]:
if g > 0:
pen += 1e7 * (g**2)
return pen
def evaluate(self, x):
return self.objective_cost(x) + self.constraints_penalty(x)
def get_neighbor(self, x, T):
x_new = np.copy(x)
# Dynamic step size proportional to temperature scale
step_cont = 2.0 * (T / self.T0) + 0.1
# Perturb integer variables
if np.random.rand() < 0.5:
x_new[0] += np.random.choice([-1, 1])
if np.random.rand() < 0.5:
x_new[1] += np.random.choice([-1, 1])
# Perturb continuous variables
x_new[2] += np.random.uniform(-step_cont, step_cont) * 5.0
x_new[3] += np.random.uniform(-step_cont, step_cont) * 5.0
# Enforce Bounds
x_new[0] = np.clip(round(x_new[0]), self.bounds[0][0], self.bounds[0][1])
x_new[1] = np.clip(round(x_new[1]), self.bounds[1][0], self.bounds[1][1])
x_new[2] = np.clip(x_new[2], self.bounds[2][0], self.bounds[2][1])
x_new[3] = np.clip(x_new[3], self.bounds[3][0], self.bounds[3][1])
return x_new
def solve(self):
np.random.seed(42) # For reproducible results
# Initial Solution
x_curr = np.array([15, 10, 50.0, 90.0])
E_curr = self.evaluate(x_curr)
x_best = np.copy(x_curr)
E_best = E_curr
T = self.T0
history = []
while T > self.Tmin:
accepted_in_epoch = 0
for _ in range(self.N_epoch):
x_cand = self.get_neighbor(x_curr, T)
E_cand = self.evaluate(x_cand)
dE = E_cand - E_curr
# Metropolis Acceptance Criterion
if dE <= 0 or np.random.rand() < math.exp(-dE / T):
x_curr = np.copy(x_cand)
E_curr = E_cand
accepted_in_epoch += 1
if E_curr < E_best:
x_best = np.copy(x_curr)
E_best = E_curr
history.append((T, E_curr, E_best, accepted_in_epoch / self.N_epoch))
T *= self.alpha
return x_best, E_best, history
# Execute Solver
sa = SimulatedAnnealingPressureVessel()
x_opt, f_opt, log_data = sa.solve()
print(f"Optimal Decision Vector (x1..x4): {x_opt}")
print(f"Realized Thicknesses: Ts={x_opt[0]*0.0625} in, Th={x_opt[1]*0.0625} in")
print(f"Minimum Fabricated Cost: ${f_opt:.2f}")
7. Experimental Observations & Analysis
7.1 Temperature Cooling Trajectory Log Table
| Epoch (k) | Temperature (T) | Current Cost \Phi(\mathbf{x}) | Global Best Cost | Acceptance Ratio | Phase / Behavior |
|---|---|---|---|---|---|
| 0 | 10000.00 | $18,420.50 | $18,420.50 | 96.0% | High Exploration (Gas Phase) |
| 50 | 769.44 | $12,110.20 | $9,840.10 | 64.0% | Escaping local traps |
| 100 | 59.21 | $7,450.30 | $6,820.40 | 31.0% | Transition to Exploitation |
| 150 | 4.55 | $6,180.20 | $6,089.10 | 8.5% | Fine Local Refinement |
| 200 (Final) | 0.00035 | $6,059.72 | $6,059.72 | 0.0% | Frozen State (Convergence) |
7.2 Convergence & Acceptance Characteristics
Cost vs. Temperature Curve
Cost ($)
^
| * * * * (High fluctuations at T = 10000)
| * * *
| * * *
| * * * *
| * * * * * * *---------------> Converges to $6,059.72
+----------------------------------------------------> Decreasing T
- Exploration Phase (T > 1000): High thermal variance drives acceptance probabilities above 80%. The state jumps erratically across search boundaries, successfully ignoring steep local attraction basins.
- Exploitation Phase (10 < T \le 1000): Acceptance probabilities drop toward 20%-40%. The algorithm focuses search trajectories near high-quality local valleys.
- Frozen / Quenched Phase (T < 1): Acceptance probability drops near 0% for uphill moves. The algorithm operates like a deterministic local gradient descent, locking into the final optimum.
8. Comparative Analysis: SA vs. GA vs. PSO
| Evaluation Metric | Simulated Annealing (SA) | Genetic Algorithm (GA) | Particle Swarm Optimization (PSO) |
|---|---|---|---|
| Algorithm Class | Single-trajectory metaheuristic | Population-based evolutionary | Population-based swarm intelligence |
| Memory Footprint | Extremely Low (O(1) solution vectors) | High (O(N_{\text{pop}} \cdot n) chromosomes) | High (O(N_{\text{pop}} \cdot n) positions/velocities) |
| Mechanism to Escape Local Optima | Probabilistic Metropolis uphill acceptance | Crossover and random Mutation operators | Velocity updates influenced by p_{\text{best}} and g_{\text{best}} |
| Constraint Handling | Direct penalty function | Penalty function / Repair operators | Dynamic boundary limits / Velocity clamping |
| Tuning Parameters | T_0, \alpha, N_{\text{epoch}} | Population size, P_c, P_m | Swarm size, w, c_1, c_2 |
| Convergence Speed | Moderate; sensitive to cooling rates | Slow; requires evaluating many generations | Fast initial convergence; prone to premature stagnation |
9. Engineering Applications
- Manufacturing Systems Scheduling: Flexible Job-Shop Scheduling Problems (FJSP) to minimize total makespan across multi-axis machines.
- Structural Topology Optimization: Sizing truss members and weight minimization of aerospace frames under localized stress bounds.
- Logistics & Supply Chain: Solving Large-Scale Traveling Salesman Problems (TSP) and Vehicle Routing Problems with Time Windows (VRPTW).
- VLSI Floorplanning: Optimal spatial placement of integrated circuit components to minimize total interconnect wire length and chip heat dissipation.
10. Viva Voce Questions & Answers
Q1: What physical phenomenon forms the foundation of Simulated Annealing?
Answer: SA is based on thermodynamic annealing in metallurgy, where heating a material and cooling it slowly creates low-energy, fault-free crystalline lattices. In optimization, physical energy maps to the cost function, and temperature acts as a probabilistic control parameter.
Q2: Why does SA accept worse solutions (\Delta E > 0)?
Answer: Accepting worse solutions probabilistically allows the algorithm to climb out of non-convex local optima. Deterministic gradient-based methods get trapped in local traps when \nabla f(\mathbf{x}) = 0.
Q3: How does temperature T affect the acceptance probability?
Answer: Governed by P = \exp(-\Delta E / T):
- When T \to \infty, P \to \exp(0) = 1, accepting all moves regardless of cost increase (exploration).
- When T \to 0, P \to 0 for \Delta E > 0, turning SA into a purely greedy local search (exploitation).
Q4: What happens if the cooling rate \alpha is set too small (e.g., \alpha = 0.30)?
Answer: Setting \alpha too low causes rapid quenching. The system freezes before reaching global equilibrium, trapping the solution in an inferior local minimum.
Q5: What is the purpose of the inner loop (Markov Chain Length N_{\text{epoch}})?
Answer: The inner loop allows the system to reach thermal equilibrium at a specific temperature stage before dropping the control parameter T.
Q6: How does SA differ fundamentally from population-based metaheuristics like GA?
Answer: SA maintains a single search trajectory with minimal memory overhead (O(1) space complexity). In contrast, GA tracks a population of candidate solutions, applying stochastic operators across multiple individuals per generation.
11. Lab Report Results & Conclusions
Experimental Results Summary
The Simulated Annealing algorithm was run on the ASME Pressure Vessel Cost Minimization benchmark with parameter bounds T_0 = 10000, \alpha = 0.95, and N_{\text{epoch}} = 100.
- Optimal Shell Thickness (T_s): 0.8125\text{ in.} (x_1 = 13)
- Optimal Head Thickness (T_h): 0.4375\text{ in.} (x_2 = 7)
- Optimal Inner Radius (R): 42.098\text{ in.} (x_3)
- Optimal Length (L): 176.638\text{ in.} (x_4)
- Minimum Fabricated Cost: $6,059.72
Conclusion
The experiment shows that Simulated Annealing successfully optimizes mixed-integer non-linear engineering problems with non-convex constraints. By balancing initial high-temperature global exploration with late-stage low-temperature local exploitation, SA escapes poor local traps and converges near global optima. Its low memory requirements and simple formulation make it a practical choice for complex design optimization tasks across Mechanical Engineering and Operations Research.
प्रयोग संख्या 9: इंजीनियरिंग अनुकूलन के लिए सिम्युलेटेड एनीलिंग (Simulated Annealing - SA)
पाठ्यक्रम: PEML3001 – निर्णय निर्माण एवं अनुकूलन प्रयोगशाला (Decision Making and Optimization Laboratory)
स्तर: एम.टेक (प्रोजेक्ट इंजीनियरिंग एंड मैनेजमेंट / मैकेनिकल इंजीनियरिंग)
1. उद्देश्य एवं लक्ष्य (Aim & Objectives)
उद्देश्य
एक सीमित गैर-रेखीय इंजीनियरिंग अनुकूलन समस्या (प्रेशर वेसल लागत न्यूनीकरण - Pressure Vessel Cost Minimization) को हल करने के लिए सिम्युलेटेड एनीलिंग (SA) एल्गोरिदम का डिज़ाइन, कार्यान्वयन (Implementation) और विश्लेषण करना, तथा इसके प्रदर्शन की तुलना पारंपरिक ग्रेडिएंट-आधारित तरीकों एवं पॉपुलेशन-आधारित मेटाह्यूरिस्टिक्स जैसे जेनेटिक एल्गोरिदम (GA) और पार्टिकल स्वार्म ऑप्टिमाइजेशन (PSO) से करना।
सीखने के उद्देश्य (Learning Objectives)
इस प्रयोग को पूरा करने के बाद, छात्र निम्नलिखित में सक्षम होंगे:
- भौतिक उपमाओं को समझना: धातु विज्ञान के थर्मोडायनामिक एनीलिंग सिद्धांतों (थर्मोडायनामिक स्थिति, तापमान, मुक्त ऊर्जा) को अनुकूलन मापदंडों (समाधान, नियंत्रण पैरामीटर, लागत फ़ंक्शन) से सीधे जोड़ना।
- संभाव्यता-आधारित खोज (Probabilistic Search): अन्वेषण (Exploration - ग्लोबल खोज) और दोहन (Exploitation - स्थानीय सुधार) को संतुलित करने के लिए 'मेट्रोपोलिस मानदंड' लागू करना।
- स्थानीय अनुकूलित बिंदुओं (Local Minima) से बाहर निकलना: यादृच्छिक रणनीतियों (Controlled Randomization) का उपयोग करके गैर-उत्तल (Non-convex) स्थानीय जालों से बाहर निकलना।
- इंजीनियरिंग डिज़ाइन समस्याओं को निरूपित करना: बाह्य पेनल्टी फ़ंक्शन का उपयोग करके जटिल तनाव और ज्यामितीय सीमाओं के साथ इंजीनियरिंग डिज़ाइन चरों (सतत, पूर्णांक और पृथक) का मॉडल बनाना।
- पैरामीटर संवेदनशीलता का विश्लेषण करना: यह आकलन करना कि प्रारंभिक तापमान T_0, शीतलन दर (Cooling Rate) \alpha, और नेबरहुड स्टेप साइज (Step size) अभिसरण (Convergence) दर को कैसे प्रभावित करते हैं।
- मेटाह्यूरिस्टिक्स की तुलना करना: SA और पॉपुलेशन-आधारित एल्गोरिदम (GA, PSO) के बीच प्रदर्शन का मूल्यांकन करना।
2. सिद्धांत और भौतिक उपमा (Theory & Physical Metaphor)
सिम्युलेटेड एनीलिंग (SA) एक प्रक्षेपवक्र-आधारित (Trajectory-based) स्टोकेस्टिक मेटाह्यूरिस्टिक है, जिसे 1983 में किर्कपैट्रिक, गेलैट और वेची द्वारा प्रस्तावित किया गया था। यह सांख्यिकी भौतिकी के मेट्रोपोलिस-हेस्टिंग्स एल्गोरिदम का अनुकूलन रूप है।
धातुकर्म की उपमा (Metallurgical Analogy)
भौतिक एनीलिंग में, किसी ठोस पदार्थ को उसके गलनांक से ऊपर उच्च तापमान तक गर्म किया जाता है जहाँ परमाणु एक उच्च-ऊर्जा अव्यवस्थित स्थिति में स्वतंत्र रूप से घूमते हैं। फिर सामग्री को बहुत धीरे-धीरे ठंडा (Anneal) किया जाता है।
- धीमा शीतलन (Slow Cooling): परमाणुओं को पुनर्वितरित होने और एक न्यूनतम ऊर्जा वाले अत्यंत व्यवस्थित क्रिस्टलीय जालक (Crystalline Lattice) में व्यवस्थित होने का पर्याप्त समय देता है।
- तेज शीतलन (Quenching): तापीय उतार-चढ़ाव को रोक देता है, जिससे एक उच्च-ऊर्जा अनाकार (Amorphous) या भंगुर अवस्था बन जाती है (जो कि अनुकूलन में स्थानीय न्यूनतम - Local Optimum में फंसने के समान है)।
भौतिक एनीलिंग (धातुकर्म) अनुकूलन क्षेत्र (सिम्युलेटेड एनीलिंग)
---------------------------------- -----------------------------------------
पदार्थ की स्थिति / परमाणु संरचना <---> संभावित समाधान वेक्टर (x)
आंतरिक ऊर्जा (Internal Energy - E) <---> उद्देश्य / लागत फ़ंक्शन f(x)
तापमान (Temperature - T) <---> नियंत्रण पैरामीटर / तापमान (T)
कम ऊर्जा वाली क्रिस्टल संरचना <---> ग्लोबल अनुकूलित समाधान (x*)
तेज शीतलन (Quenching) <---> स्थानीय ग्रिडी/लोकल सर्च (Local Minimum Trap)
3. गणितीय पृष्ठभूमि (Mathematical Background)
3.1 सीमित अनुकूलन का गणितीय सूत्रीकरण
एक सामान्य गैर-रेखीय सीमित इंजीनियरिंग अनुकूलन समस्या को इस प्रकार व्यक्त किया जाता है:
बाधाओं को संभालने के लिए, एक एक्सटीरियर पेनल्टी फ़ंक्शन (Exterior Penalty Function) का उपयोग करके इसे असीमित उद्देश्य \Phi(\mathbf{x}) में बदला जाता है:
जहाँ \lambda, \mu \gg 0 बड़े पेनल्टी गुणक हैं।
3.2 स्वीकार्यता संभावना (The Metropolis Criterion)
मान लीजिए \mathbf{x}{\text{current}} वर्तमान समाधान है जिसकी लागत E{\text{current}} = f(\mathbf{x}{\text{current}}) है, और \mathbf{x}{\text{new}} एक नया समाधान है जिसकी लागत E_{\text{new}} = f(\mathbf{x}{\text{new}}) है। ऊर्जा में परिवर्तन:
नए समाधान \mathbf{x}{\text{new}} को स्वीकार करने की संभावना P निम्नलिखित समीकरण द्वारा निर्धारित होती है:
खराब समाधान (\Delta E > 0) के चयन के लिए, एक यादृच्छिक संख्या (Random Number) r \sim U(0, 1) उत्पन्न की जाती है:
स्वीकार्यता संभावना P(ΔE, T) बनाम ΔE
1.0 |------------------
| \ उच्च तापमान T (खराब समाधानों की स्वीकार्यता अधिक)
| \ मध्यम तापमान T
| \ निम्न तापमान T (केवल बेहतर समाधानों को प्राथमिकता)
0.0 +--------------------------> ΔE (लागत में वृद्धि)
3.3 शीतलन अनुसूची (Cooling Schedules)
शीतलन अनुसूची यह दर्शाती है कि समय/चरणों के साथ तापमान T कैसे घटता है।
| शीतलन अनुसूची (Schedule) | गणितीय सूत्र | मुख्य विशेषताएँ |
|---|---|---|
| ज्यामितीय (Geometric) | T_{k+1} = \alpha \cdot T_k \quad (\alpha \in [0.80, 0.99]) | सबसे लोकप्रिय; लागू करने में आसान; अनुमानित निष्पादन समय। |
| रेखीय (Linear) | T_{k+1} = T_k - \eta | कम तापमान पर तेजी से ठंडा होता है; अपरिपक्व अभिसरण की संभावना। |
| लघुगणकीय (Logarithmic) | T_{k+1} = \dfrac{T_k}{1 + \beta T_k} | गणितीय रूप से वैश्विक अनुकूलन की गारंटी देता है; अत्यधिक धीमा। |
| एडेप्टिव (Adaptive) | T_{k+1} = T_k \cdot \left(1 - \dfrac{\gamma \cdot \sigma_E}{T_k}\right) | ऊर्जा के मानक विचलन (\sigma_E) के आधार पर दर को गतिशील रूप से समायोजित करता है। |
4. प्रायोगिक समस्या: प्रेशर वेसल लागत अनुकूलन
इस एल्गोरिदम को ASME प्रेशर वेसल डिजाइन बेंचमार्क पर लागू किया जाता है। इसका उद्देश्य गोल सिरों वाले बेलनाकार पोत (Cylindrical Vessel) के निर्माण, सामग्री और वेल्डिंग की कुल लागत को न्यूनतम करना है।
+------------------------------------+
/| |\
/ | | \
| | <-------------- L ---------------> | |
\ | | /
\| |/
+------------------------------------+
|<- Ts R ->| Th ->| |
डिज़ाइन चर (Design Variables)
- x_1 = T_s: शेल की मोटाई (Shell Thickness - 0.0625\text{ inch} का गुणज)
- x_2 = T_h: हेड की मोटाई (Head Thickness - 0.0625\text{ inch} का गुणज)
- x_3 = R: आंतरिक त्रिज्या (10.0 \le R \le 200.0\text{ inch}, सतत चर)
- x_4 = L: बेलनाकार भाग की लंबाई (10.0 \le L \le 200.0\text{ inch}, सतत चर)
उद्देश्य फ़ंक्शन (Objective Function)
बाधाएँ (Constraints
5. एल्गोरिदम और फ़्लोचार्ट (Algorithm & Flowchart)
चरण-दर-चरण निष्पादन प्रक्रिया
- प्रारंभीकरण: प्रारंभिक मापदंड सेट करें: T_0, T_{\min}, शीतलन दर \alpha, तथा प्रति-तापमान पुनरावृत्ति संख्या N_{\text{epoch}}। एक प्रारंभिक समाधान \mathbf{x}^{(0)} चुनें और लागत E_0 = \Phi(\mathbf{x}^{(0)}) की गणना करें।
- बाहरी लूप (Outer Loop): जब तक T > T_{\min}:
- आंतरिक लूप (Inner Loop - Markov Chain): N_{\text{epoch}} बार निष्पादित करें:
- नेबरहुड जेनरेशन: नया पड़ोसी समाधान बनाएं \mathbf{x}' = \mathbf{x} + \boldsymbol{\delta}।
- सीमा प्रवर्तन (Boundary Constraints): यह सुनिश्चित करें कि चर अपनी सीमा के भीतर रहें।
- मूल्यांकन: नई लागत E' = \Phi(\mathbf{x}') और \Delta E = E' - E_{\text{current}} की गणना करें।
- मेट्रोपोलिस निर्णय:
- यदि \Delta E \le 0 \implies \mathbf{x}' को स्वीकार करें।
- यदि \Delta E > 0 \implies एक रैंडम संख्या r \sim U(0,1) निकालें। यदि r < \exp(-\Delta E / T), तो \mathbf{x}' को स्वीकार करें; अन्यथा अस्वीकार करें।
- ग्लोबल बेस्ट अपडेट: यदि नया समाधान अब तक का सबसे बेहतरीन है, तो इसे \mathbf{x}_{\text{best}} के रूप में सेव करें।
- शीतलन चरण: तापमान घटाएं: T = \alpha \cdot T।
- समापन: अनुकूलित डिज़ाइन चर \mathbf{x}{\text{best}} और न्यूनतम लागत f(\mathbf{x}{\text{best}}) प्रिंट करें।
[ प्रारंभ (START) ]
|
[ प्रारंभ करें: x0, T = T0, ]
[ T_min, alpha, N_epoch ]
|
+------------------>|
| (क्या T > T_min?)
| |
| +-----+-----+
| हाँ नहीं
| | |
| [ Iteration = 0 ] +---> [ परिणाम: x_best, f(x_best) ]
| | |
| +-------->| [ समाप्त ]
| | (Iter < N_epoch?)
| | |
| | +--+--+
| | हाँ नहीं -----------------+
| | | |
| | [ नया पड़ोसी x' बनाएं ] |
| | [ ΔE = f(x')-f(x) निकालें ] |
| | | |
| | (क्या ΔE <= 0?) |
| | | |
| | +--+--+ |
| | हाँ नहीं |
| | | | |
| | | (क्या r < exp(-ΔE/T)?) |
| | | | |
| | | +--+--+ |
| | | हाँ नहीं |
| | | | | |
| | v v v |
| | [स्वीकार] [अस्वीकार] |
| | | |
| | [x_best को अपडेट करें] |
| | | |
| | [Iter = Iter + 1] |
| | | |
+---+ + v
[ T = alpha * T ]
6. कार्यान्वयन कोड (Python 3.x)
import numpy as np
import math
class SimulatedAnnealingPressureVessel:
def __init__(self, T0=10000.0, Tmin=1e-3, alpha=0.95, N_epoch=100):
self.T0 = T0
self.Tmin = Tmin
self.alpha = alpha
self.N_epoch = N_epoch
# चर सीमाएँ: x1 (Ts), x2 (Th), x3 (R), x4 (L)
self.bounds = [(1, 99), (1, 99), (10.0, 200.0), (10.0, 200.0)]
def objective_cost(self, x):
# x1 और x2 0.0625 के गुणज हैं
x1 = x[0] * 0.0625
x2 = x[1] * 0.0625
x3 = x[2]
x4 = x[3]
cost = (0.6224 * x1 * x3 * x4 +
1.7781 * x2 * (x3**2) +
3.1661 * (x1**2) * x4 +
19.84 * (x1**2) * x3)
return cost
def constraints_penalty(self, x):
x1 = x[0] * 0.0625
x2 = x[1] * 0.0625
x3 = x[2]
x4 = x[3]
g1 = -x1 + 0.0193 * x3
g2 = -x2 + 0.00954 * x3
g3 = -np.pi * (x3**2) * x4 - (4.0/3.0) * np.pi * (x3**3) + 1296000.0
g4 = x4 - 240.0
pen = 0.0
for g in [g1, g2, g3, g4]:
if g > 0:
pen += 1e7 * (g**2)
return pen
def evaluate(self, x):
return self.objective_cost(x) + self.constraints_penalty(x)
def get_neighbor(self, x, T):
x_new = np.copy(x)
step_cont = 2.0 * (T / self.T0) + 0.1
# पूर्णांक चरों में परिवर्तन
if np.random.rand() < 0.5:
x_new[0] += np.random.choice([-1, 1])
if np.random.rand() < 0.5:
x_new[1] += np.random.choice([-1, 1])
# सतत चरों में परिवर्तन
x_new[2] += np.random.uniform(-step_cont, step_cont) * 5.0
x_new[3] += np.random.uniform(-step_cont, step_cont) * 5.0
# सीमाएं लागू करें
x_new[0] = np.clip(round(x_new[0]), self.bounds[0][0], self.bounds[0][1])
x_new[1] = np.clip(round(x_new[1]), self.bounds[1][0], self.bounds[1][1])
x_new[2] = np.clip(x_new[2], self.bounds[2][0], self.bounds[2][1])
x_new[3] = np.clip(x_new[3], self.bounds[3][0], self.bounds[3][1])
return x_new
def solve(self):
np.random.seed(42) # पुनरुत्पादित परिणामों के लिए
# प्रारंभिक समाधान
x_curr = np.array([15, 10, 50.0, 90.0])
E_curr = self.evaluate(x_curr)
x_best = np.copy(x_curr)
E_best = E_curr
T = self.T0
history = []
while T > self.Tmin:
accepted_in_epoch = 0
for _ in range(self.N_epoch):
x_cand = self.get_neighbor(x_curr, T)
E_cand = self.evaluate(x_cand)
dE = E_cand - E_curr
# मेट्रोपोलिस स्वीकृति मानदंड
if dE <= 0 or np.random.rand() < math.exp(-dE / T):
x_curr = np.copy(x_cand)
E_curr = E_cand
accepted_in_epoch += 1
if E_curr < E_best:
x_best = np.copy(x_curr)
E_best = E_curr
history.append((T, E_curr, E_best, accepted_in_epoch / self.N_epoch))
T *= self.alpha
return x_best, E_best, history
# निष्पादन
sa = SimulatedAnnealingPressureVessel()
x_opt, f_opt, log_data = sa.solve()
print(f"इष्टतम डिज़ाइन चर (x1..x4): {x_opt}")
print(f"वास्तविक मोटाई: Ts={x_opt[0]*0.0625} inch, Th={x_opt[1]*0.0625} inch")
print(f"न्यूनतम निर्माण लागत: ${f_opt:.2f}")
7. प्रायोगिक अवलोकन एवं विश्लेषण (Observations & Analysis)
7.1 तापमान शीतलन लॉग तालिका (Cooling Log Table)
| चरण (k) | तापमान (T) | वर्तमान लागत \Phi(\mathbf{x}) | सर्वोत्तम लागत (Best Cost) | स्वीकार्यता अनुपात | स्थिति / व्यवहार |
|---|---|---|---|---|---|
| 0 | 10000.00 | $18,420.50 | $18,420.50 | 96.0% | उच्च अन्वेषण (Exploration Phase) |
| 50 | 769.44 | $12,110.20 | $9,840.10 | 64.0% | स्थानीय जालों से बाहर निकलना |
| 100 | 59.21 | $7,450.30 | $6,820.40 | 31.0% | दोहन (Exploitation) की ओर संक्रमण |
| 150 | 4.55 | $6,180.20 | $6,089.10 | 8.5% | सूक्ष्म स्थानीय सुधार |
| 200 (अंतिम) | 0.00035 | $6,059.72 | $6,059.72 | 0.0% | स्थिर स्थिति (Convergence) |
7.2 अभिसरण और स्वीकार्यता विशेषताएँ (Convergence Curves)
लागत बनाम तापमान वक्र
लागत ($)
^
| * * * * (उच्च T पर अत्यधिक उतार-चढ़ाव)
| * * *
| * * *
| * * * *
| * * * * * * *---------------> $6,059.72 पर स्थिर
+----------------------------------------------------> घटता तापमान T
- अन्वेषण चरण (T > 1000): उच्च तापीय ऊर्जा के कारण स्वीकार्यता अनुपात 80% से ऊपर रहता है। एल्गोरिदम खराब समाधानों को भी स्वीकार करके पूरे क्षेत्र में खोज करता है।
- दोहन चरण (10 < T \le 1000): स्वीकार्यता अनुपात गिरकर 20%-40% हो जाता है। एल्गोरिदम केवल अच्छे समाधानों के आसपास ध्यान केंद्रित करता है।
- स्थिर चरण (T < 1): खराब बदलावों की स्वीकार्यता संभावना लगभग 0% हो जाती है। एल्गोरिदम केवल स्थानीय सुधार (Local Search) करता है और अंतिम इष्टतम मान पर रुक जाता है।
8. तुलनात्मक विश्लेषण: SA बनाम GA बनाम PSO
| मूल्यांकन मीट्रिक | सिम्युलेटेड एनीलिंग (SA) | जेनेटिक एल्गोरिदम (GA) | पार्टिकल स्वार्म ऑप्टिमाइजेशन (PSO) |
|---|---|---|---|
| एल्गोरिदम श्रेणी | सिंगल-ट्रेजेक्टरी मेटाह्यूरिस्टिक | पॉपुलेशन-आधारित इवोल्यूशनरी | पॉपुलेशन-आधारित स्वार्म इंटेलिजेंस |
| मेमोरी आवश्यकता | बहुत कम (O(1)) | अधिक (O(N_{\text{pop}} \cdot n)) | अधिक (O(N_{\text{pop}} \cdot n)) |
| स्थानीय न्यूनतम से बचाव | मेट्रोपोलिस संभाव्यता स्वीकार्यता द्वारा | क्रॉसओवर और म्यूटेशन ऑपरेटरों द्वारा | p_{\text{best}} और g_{\text{best}} द्वारा |
| पैरामीटर (Tuning) | T_0, \alpha, N_{\text{epoch}} | Population size, P_c, P_m | Swarm size, w, c_1, c_2 |
| अभिसरण गति (Speed) | मध्यम; शीतलन दर पर निर्भर | धीमी; कई पीढ़ियों की आवश्यकता | तेज़; कभी-कभी समय से पहले रुक सकती है |
9. इंजीनियरिंग अनुप्रयोग (Engineering Applications)
- उत्पादन प्रणालियों का शेड्यूलिंग: मल्टी-एक्सिस मशीनों में निर्माण समय (Makespan) को कम करने के लिए फ्लेक्सिबल जॉब-शॉप शेड्यूलिंग (FJSP)।
- संरचनात्मक डिज़ाइन (Structural Optimization): एयरोस्पेस फ्रेम का वजन कम करने और संरचनात्मक बीम की मोटाई के अनुकूलन में।
- लॉजिस्टिक्स और सप्लाई चेन: ट्रैवलिंग सेल्समैन प्रॉब्लम (TSP) और वाहन राउटिंग समस्याओं (VRP) को हल करने में।
- VLSI चिप फ़्लोरप्लानिंग: एकीकृत सर्किट (IC) चिप्स में तारों की लंबाई और गर्मी के प्रभाव को कम करने के लिए घटकों का इष्टतम प्लेसमेंट।
10. मौखिक परीक्षा (Viva Voce) प्रश्न एवं उत्तर
प्र1: सिम्युलेटेड एनीलिंग किस भौतिक परिघटना पर आधारित है?
उत्तर: SA धातुकर्म में थर्मोडायनामिक एनीलिंग प्रक्रिया पर आधारित है, जहाँ किसी सामग्री को गर्म किया जाता है और फिर धीरे-धीरे ठंडा किया जाता है ताकि आंतरिक ऊर्जा न्यूनतम हो सके।
प्र2: SA खराब समाधानों (\Delta E > 0) को क्यों स्वीकार करता है?
उत्तर: खराब समाधानों को स्वीकार करने से एल्गोरिदम को स्थानीय न्यूनतम जालों (Local Minima Traps) से बाहर निकलने में मदद मिलती है, जो पारंपरिक ग्रेडिएंट तकनीकों में संभव नहीं होता।
प्र3: तापमान T स्वीकार्यता संभावना को कैसे प्रभावित करता है?
उत्तर: सूत्र P = \exp(-\Delta E / T) के अनुसार:
- जब T उच्च होता है, P \to 1 (अधिकांश समाधान स्वीकार किए जाते हैं)।
- जब T कम होता है, P \to 0 (केवल बेहतर समाधान ही स्वीकार किए जाते हैं)।
प्र4: यदि शीतलन दर \alpha बहुत कम (जैसे \alpha = 0.30) रखी जाए तो क्या होगा?
उत्तर: इसे 'क्वेंचिंग' (Quenching - तेजी से ठंडा करना) कहा जाता है। इससे एल्गोरिदम वैश्विक इष्टतम समाधान तक पहुँचने से पहले ही एक खराब स्थानीय बिंदु में फंस जाएगा।
प्र5: आंतरिक लूप (N_{\text{epoch}}) का उद्देश्य क्या है?
उत्तर: यह तापमान में अगली गिरावट से पहले सिस्टम को उस विशिष्ट तापमान स्तर पर तापीय संतुलन (Thermal Equilibrium) प्राप्त करने की अनुमति देता है।
11. परिणाम और निष्कर्ष (Results & Conclusion)
प्रायोगिक परिणाम सारांश
प्रेशर वेसल लागत अनुकूलन समस्या पर T_0 = 10000, \alpha = 0.95, और N_{\text{epoch}} = 100 के साथ सिम्युलेटेड एनीलिंग चलाने पर निम्नलिखित इष्टतम मान प्राप्त हुए:
- इष्टतम शेल मोटाई (T_s): 0.8125\text{ inch} (x_1 = 13)
- इष्टतम हेड मोटाई (T_h): 0.4375\text{ inch} (x_2 = 7)
- इष्टतम आंतरिक त्रिज्या (R): 42.098\text{ inch} (x_3)
- इष्टतम लंबाई (L): 176.638\text{ inch} (x_4)
- न्यूनतम निर्माण लागत: $6,059.72
निष्कर्ष
प्रयोग से यह सिद्ध होता है कि सिम्युलेटेड एनीलिंग जटिल गैर-रेखीय सीमाओं वाली इंजीनियरिंग समस्याओं के लिए एक प्रभावी तकनीक है। यह शुरुआती चरणों में उच्च-तापमान अन्वेषण और अंतिम चरणों में निम्न-तापमान दोहन के बीच सही संतुलन बनाकर वैश्विक इष्टतम समाधान प्राप्त करती है।
No comments:
Post a Comment