Weld Beam Design Optimization¶

"""
Welded Beam Design Problem
Example of constrained optimization using GA with PyEGRO

Objective:
    Minimize fabrication cost of a welded beam
Constraints:
    - Shear stress ≤ 13,600 psi
    - Normal stress ≤ 30,000 psi
    - Buckling, deflection, geometry, and minimum thickness
"""

import numpy as np
from PyEGRO.optimize.GA import run_deterministic_optimization, save_optimization_results

# === Problem constants ===
P = 6000         # load in lb
L = 14           # length in inches
E = 30e6         # psi
G = 12e6         # psi
tau_max = 13600
sigma_max = 30000
delta_max = 0.25

# === Objective function ===
def weld_beam_cost(X):
    X = np.atleast_2d(X)
    return np.array([
        1.10471 * h**2 * l + 0.04811 * t * b * (14 + l)
        for h, l, t, b in X
    ])

# === Constraint functions ===
def constraint_tau(X):
    results = []
    for h, l, t, b in np.atleast_2d(X):
        R = np.sqrt(l**2 / 4 + ((h + t) / 2)**2)
        J = 2 * (np.sqrt(2) * h * l * ((l**2 / 12) + ((h + t) / 2)**2))
        tau_p = P / (np.sqrt(2) * h * l)
        tau_d = (P * (L + l / 2) * R) / J
        tau = np.sqrt(tau_p**2 + 2 * tau_p * tau_d * l / (2 * R) + tau_d**2)
        results.append(tau - tau_max)
    return np.array(results)

def constraint_sigma(X):
    return np.array([
        6 * P * L / (b * t**2) - sigma_max
        for h, l, t, b in np.atleast_2d(X)
    ])

def constraint_geometry(X):
    return np.array([h - b for h, l, t, b in np.atleast_2d(X)])  # h - b <= 0

def constraint_buckling(X):
    results = []
    for h, l, t, b in np.atleast_2d(X):
        Pc = (4.013 * E * np.sqrt(t**2 * b**6 / 36)) / (L**2) * (1 - t / (2 * L) * np.sqrt(E / (4 * G)))
        results.append(P - Pc)
    return np.array(results)

def constraint_deflection(X):
    return np.array([
        4 * P * L**3 / (E * b * t**3) - delta_max
        for h, l, t, b in np.atleast_2d(X)
    ])

def constraint_thickness(X):
    return np.array([0.125 - h for h, l, t, b in np.atleast_2d(X)])

# === Define problem info ===
data_info = {
    'variables': [
        {'name': 'h', 'vars_type': 'design_vars', 'range_bounds': [0.125, 5], 'description': 'Weld thickness'},
        {'name': 'l', 'vars_type': 'design_vars', 'range_bounds': [0.1, 10], 'description': 'Length of welded joint'},
        {'name': 't', 'vars_type': 'design_vars', 'range_bounds': [0.1, 10], 'description': 'Width of the beam'},
        {'name': 'b', 'vars_type': 'design_vars', 'range_bounds': [0.1, 5], 'description': 'Thickness of the beam'}
    ]
}

# === List of constraint functions ===
constraint_functions = [
    constraint_tau,
    constraint_sigma,
    constraint_geometry,
    constraint_buckling,
    constraint_deflection,
    constraint_thickness
]

# === Run optimization ===
results = run_deterministic_optimization(
    data_info=data_info,
    true_func=weld_beam_cost,
    constraint_funcs=constraint_functions,
    pop_size=600,
    n_gen=100,
    sampling_method='lhs',
    crossover_prob=0.9,
    crossover_eta=15,
    mutation_eta=20
)

# === Save results ===
save_optimization_results(
    results=results,
    data_info=data_info,
    save_dir='WELD_BEAM_DIRECT_CONSTRAINED'
)

# === Print result ===
print("\nOptimized Solution:")
h, l, t, b = results['best_solution']
print(f"  h: {h:.6f}")
print(f"  l: {l:.6f}")
print(f"  t: {t:.6f}")
print(f"  b: {b:.6f}")
print(f"Objective Value: {results['best_fitness']:.6f}")
print(f"Feasible: {results['is_feasible']}")