The xlwings Scipy spreadsheet has been updated with a new example of the xl_SolveF function, that uses the Scipy Optimize root function. The new spreadsheet can be downloaded from:
The new example uses a python function ic_calc (included in the download file) to find the ultimate reaction forces and moment for a group of bolts with specified load eccentricity and angle, and yield displacement and force. The calculation procedure is detailed at:
The shear force on each bolt is non-linear, depending on the bolt strength and shear displacement, and the displacement in each bolt depends on the position of the centre of rotation. The Scipy root function adjusts the centre of rotation coordinates and the maximum bolt displacement so that the reaction forces and moments are equal to the applied forces and moment:
Results can be compared with a VBA spreadsheet doing the same calculation but using the Excel Solver, from Yakpol’s Spreadsheet Solutions for Structural Engineering.
Note that to generate the same result as my spreadsheet
- The XY coordinates must be adjusted to give the same perpendicular eccentricity from the centroid of the bolts to the line of action of the force.
- The force angle is the angle to the x-axis, whereas my function uses the angle to the Y axis.
The xl_SolveF function passes data to the function to solved as two arguments:
- A 1D vector containing the variables
- A single array containing all other data
An interface function is therefore used to extract the input parameters for the ic_calc function in the correct format, as shown below:
def ic_check(IC, vals): xloc = vals yloc = vals vals2 =vals res1 = xlic_calc(IC,xloc, yloc, vals2) totX = res1 totY = res1 totM = res1 Pu =IC Pux = res1 Puy = res1 Mu = res1 return [totX-Pux, totY-Puy, totM-Mu] @xw.func @xw.ret(transpose = True) def xlic_calc(IC, xloc, yloc, vals): xloc = np.array(xloc) yloc = np.array(yloc) ecc =vals theta = np.radians(vals) deltamax = vals Rult = vals num_bolts = len(xloc) ICx = IC ICy = IC Pu = IC Pux = Pu * np.sin(theta) Puy = Pu * np.cos(theta) Mu = Pu * (ecc-ICx*np.cos(theta)-ICy*np.sin(theta)) xIC = xloc - ICx yIC = yloc - ICy di = np.sqrt((xIC*xIC)+(yIC*yIC)) dmax = max(di) deltai = di/dmax * deltamax ri = Rult * (1-np.e**(-10.0*deltai))**0.55 fx = yIC * ri / di fy = xIC * ri / di moment = ri * di totX = sum(fx) totY = sum(fy) totM = sum(moment) return [totX, totY, totM, Pux, Puy, Mu]