Optimization example using PyMAPDL and SciPy Optimize library

Below is a PyMAPDL example using FLUID218 elements (3D Hydrodynamic Bearing Element) on how you could find an optimized Inputs (Shaft center X and Y location) for a preferred FLUID forces in X and Y directions.

Please refer FLUID218 Documentation below
https://ansyshelp.ansys.com/account/secured?returnurl=/Views/Secured/corp/v241/en/ans_elem/Hlp_E_FLUID218.html?q=FLUID218

Find more posts tagged with

PyMAPDL

MAPDL

Comments

Vishnu

The first step in optimization is formulating an objective function, in this case the function named extract_fluid_forces is evaluated by using PyMAPDL

scipy.optimize ,minimize is used for optimization .please refer the documentation for it below:

https://docs.scipy.org/doc/scipy/reference/generated/scipy.optimize.minimize.html

Example code below

# Import necessary libraries
import matplotlib.pyplot as plt
import numpy as np
from ansys.dpf import core as dpf
from scipy.optimize import minimize

from ansys.mapdl.core import launch_mapdl

# Define the path to the result file
rstpath = r"D:\apdl\test\file.rst"

# Launch MAPDL with specified run location and number of processors
mapdl = launch_mapdl(run_location=r"D:\apdl\test",nproc=1,override=True)

mapdl.clear("All")
mapdl.prep7()
mapdl.et('', 'FLUID218')
mapdl.mp(lab='VISC', mat='1', c0='0.135')
mapdl.mp(lab='DENS', mat='1', c0='910')
mapdl.r(nset='1', r1='0', r2='10')
mapdl.cyl4(xcenter='', ycenter='', rad1='10', theta1='0', rad2='', theta2='90', depth='50')
mapdl.cyl4(xcenter='', ycenter='', rad1='10', theta1='90', rad2='', theta2='180', depth='50')
mapdl.cyl4(xcenter='', ycenter='', rad1='10', theta1='180', rad2='', theta2='270', depth='50')
mapdl.cyl4(xcenter='', ycenter='', rad1='10', theta1='270', rad2='', theta2='360', depth='50')
mapdl.vdele(nv1='ALL', nv2='', ninc='', kswp='0')
mapdl.run('CSYS,1')
mapdl.asel('U', 'LOC', 'X', 9, 10)
mapdl.adele('ALL')
mapdl.allsel()
mapdl.mshape(0, '2D')
mapdl.mshkey(1)
mapdl.lesize('ALL', size='', angsiz='',  ndiv=20)
mapdl.amesh('ALL')
mapdl.eplot()
mapdl.nummrg('ALL')
clearance = 0.1
mapdl.rmodif(1, 8, clearance)
mapdl.rmodif(1, 9, clearance)
mapdl.rmodif(1, 10, clearance)
mapdl.rmodif(1, 11, clearance)
mapdl.nsel('S', 'LOC', 'Z', '', 0)
mapdl.d('ALL', 'PRES', 0)
mapdl.allsel()
mapdl.nsel('S', 'LOC', 'Z', 50, 50)
mapdl.d('ALL', 'PRES', 100)
mapdl.nsel('S', 'LOC', 'Y', -20, 20)
mapdl.d('ALL', 'PRES', 100)
mapdl.nsel('S', 'LOC', 'Y', 160, 200)
mapdl.d('ALL', 'PRES', 100)
mapdl.allsel()
mapdl.omega('', '', 10)
mapdl.finish()

# Define a function to extract fluid forces
def extract_fluid_forces(xoffset,yoffset):
        mapdl.prep7()
        mapdl.rmodif(1, 3,xoffset)
        mapdl.rmodif(1, 4,yoffset)
        mapdl.solution()
        mapdl.solve()       
        # data source
        my_data_sources = dpf.DataSources(rstpath)

        # model
        my_model = dpf.Model(rstpath)

        # mesh
        my_mesh = my_model.metadata.meshed_region

        # Field Containers at all mesh
        # Get result operator
        result_operator = dpf.operators.result.__getattribute__("nmisc")

        result_field_container = result_operator(
            time_scoping=[1], data_sources=my_data_sources, mesh=my_mesh,item_index=11
        ).outputs.fields_container.get_data()

        # First Field usually field corresponding to first time scoping
        resultfield = result_field_container[0]

        FXTOTAL=sum(resultfield.data_as_list)

        result_field_container = result_operator(
            time_scoping=[1], data_sources=my_data_sources, mesh=my_mesh,item_index=12
        ).outputs.fields_container.get_data()

        # First Field usually field corresponding to first time scoping
        resultfield = result_field_container[0]

        FYTOTAL=sum(resultfield.data_as_list)

        fluid_forces = np.array([FXTOTAL,FYTOTAL])
        #plt.bar(range(len(fluid_forces)), fluid_forces)

        return fluid_forces


# preferred forces
preferred_fx = 10000  # replace with your preferred FXtotal
preferred_fy = 90000  # replace with your preferred FYtotal

# Define objective function
def objective(offsets):
    try:
        [fx, fy] = extract_fluid_forces(offsets[0],offsets[1])
        return (((fx - preferred_fx)**2 + (fy - preferred_fy)**2))
    except:  # catch all exceptions
        return 9E50  # return a large value


# initial guess for x and y offsets
x0 = [0.0,0.0]

# Minimize the objective function
result = minimize(objective, x0,method='nelder-mead',options={'maxfev':1000,'maxiter':1000,'disp': True,'fatol' : 1e-10, 'xatol' : 1e-4})

# Get the optimized x and y offsets
x, y = result.x
fx,fy = extract_fluid_forces(x, y)

# Exit MAPDL
mapdl.exit()

# Print the optimized values
print("Optimized X and Y Offsets are %f and %f"%(x,y))
print("The Fx,Fy forces at these X and Y offsets are %f and %f "%(fx,fy))


# Define the labels for the bars
labels = ['FX', 'FY']

# Create a figure and a set of subplots
fig, ax = plt.subplots()

# Define the width of the bars
width = 0.35

# Define the positions for the groups of bars
x = np.arange(len(labels))

# Plot the preferred and obtained values as bars next to each other
ax.bar(x - width/2, [preferred_fx,preferred_fy], width, label='Preferred')
ax.bar(x + width/2, [fx,fy], width, label='optimized')

# Add some text for labels, title and custom x-axis tick labels, etc.
ax.set_ylabel('Values')
ax.set_title('Preferred and Obtained Values of FX and FY')
ax.set_xticks(x)
ax.set_xticklabels(labels)
ax.legend()


# Display the plot
plt.show()

Vishnu

Output --->>>