1

I am trying to use Gekko in optimal design of a reinforced concrete building modeled with Abaqus.

*my design variables are the cross sections of the the building's structural system. Every variable of this kind is a set of:

  • dimensions of section
  • Reinforcement of the section
  • Materials properties of the section
  • the cost of unit length of the section

*the objective function which is the total cost of the building is calculated by multiplying every section's unit length's cost by its length then summing up all these multiplies.

the constraints of the problem are calculated after running an analysis on Abaqus and then obtaining all the required results that are needed to fulfil design requirements. I mean that any break for any requirement (eg. maximum stresses or strains or deformations....) will initiate a penalty that would be summed up with the objective function and increase the cost of building significantly.

My main question is that after I created the database for my building's sections consisting of all possible sections, can I make Gekko pick a value of design variable from this database and from it only? knowing that every variable is a set of many parameters (dimensions, material props, etc....)

Christian
  • 4,902
  • 4
  • 24
  • 42

1 Answers1

0

Design optimization is possible with gekko but the equations are needed for optimization. A possible way to interface pre-computed solutions is to use csplines (1 design variable), bsplines (2 design variables), or machine learning (3 variables+). If it needs to select from a set of options then an sos1 (special-ordered set) may be needed.

Here is one example optimization problem with a tubular column.

tubular column

from gekko import GEKKO
import matplotlib.pyplot as plt
import numpy as np
from scipy.optimize import fsolve

m = GEKKO()

#%% Constants
pi = m.Const(3.14159,'pi')
P = 2300 # compressive load (kg_f)
o_y = 450 # allowable yield stress (kg_f/cm^2)
E = 0.65e6 # elasticity (kg_f/cm^2)
p = 0.0020 # weight density (kg_f/cm^3)
l = 300 # length of the column (cm)

#%% Variables (the design variables available to the solver)
d = m.Var(value=8.0,lb=2.0,ub=14.0) # mean diameter (cm)
t = m.Var(value=0.3,lb=0.2 ,ub=0.8) # thickness (cm)
cost = m.Var()

#%% Intermediates (computed by solver from design variables and constants)
d_i = m.Intermediate(d - t)
d_o = m.Intermediate(d + t)
W = m.Intermediate(p*l*pi*(d_o**2 - d_i**2)/4) # weight (kgf)
o_i = m.Intermediate(P/(pi*d*t)) # induced stress

# second moment of area of the cross section of the column
I = m.Intermediate((pi/64)*(d_o**4 - d_i**4))

# buckling stress (Euler buckling load/cross-sectional area)
o_b = m.Intermediate((pi**2*E*I/l**2)*(1/(pi*d*t)))

#%% Equations (constraints, etc. Cost could be an intermediate variable)
m.Equations([
o_i - o_y <= 0,
o_i - o_b <= 0,
cost == 5*W + 2*d
])

#%% Objective
m.Minimize(cost)

#%% Solve and print solution
m.options.SOLVER = 1
m.solve()

print('Optimal cost: ' + str(cost[0]))
print('Optimal mean diameter: ' + str(d[0]))
print('Optimal thickness: ' + str(t[0]))

minima = np.array([d[0], t[0]]) 

#%% Contour plot
# create a cost function as a function of the design variables d and t
f = lambda d, t: 2 * d + 5 * p * l * np.pi * ((d+t)**2 - (d-t)**2)/4 

xmin, xmax, xstep = 2, 14, .2 # diameter
ymin, ymax, ystep = .2, .8, .05 # thickness

d, t = np.meshgrid(np.arange(xmin, xmax + xstep, xstep), \
                   np.arange(ymin, ymax + ystep, ystep))
z = f(d, t)

# Determine the compressive stress constraint line. 
#stress = P/(pi*d*t) # induced axial stress
t_stress = np.arange(ymin, ymax, .025) # use finer step to get smoother constraint line
d_stress = []
for tt in t_stress:
    dd = P/(np.pi * tt * o_y)
    d_stress.append(dd)

# Determine buckling constraint line. This is tougher because we cannot
#  solve directly for t from d. Used scipy.optimize.fsolve to find roots 
d_buck = []
t_buck = []
for d3 in np.arange(6, xmax, .005): 
    fb = lambda t : o_y-np.pi**2*E*((d3+t)**4-(d3-t)**4)/(64*l**2*d3*t)
    tr = np.array([0.3])
    roots = fsolve(fb, tr)
    if roots[0] != 0: 
        if roots[0] >= .1 and roots[0]<=1.:
            t_buck.append(roots[0])
            d_buck.append(d3)

# Create contour plot
plt.style.use('ggplot') # to make prettier plots
fig, ax = plt.subplots(figsize=(10, 6))

CS = ax.contour(d, t, z, levels=15,) 
ax.clabel(CS, inline=1, fontsize=10)
ax.set_xlabel('mean diameter $d$')
ax.set_ylabel('half thickness $t$')

ax.set_xlim((xmin, xmax))
ax.set_ylim((ymin, ymax))

# Add constraint lines and optimal marker
ax.plot(d_stress, t_stress, "->", label="Stress constraint")
ax.plot(d_buck, t_buck, "->", label="Buckling constraint" )

minima_ = minima.reshape(-1, 1)
ax.plot(*minima_, 'r*', markersize=18, label="Optimum")
ax.text(10,.25,"Contours = Cost (objective)\nConstraint line markers point\ntowards feasible space.")
plt.title('Column Design')
plt.legend()
plt.show()
TexasEngineer
  • 684
  • 6
  • 15