mdd/lib/test/server/gfk_plate.py

# https://akaszynski.github.io/pyansys/ansys_control.html#initial-setup-and-example
# pip install pyansys --upgrade
import pyansys
import numpy as np
import os
import time
from scipy import sparse
from scipy.sparse import linalg


def sim_rotor(materials, angles, outer_diameter=0.250, inner_diameter=0.050, thickness=0.005, omega=1):
    layers = len(angles)
    layer_thickness = thickness / layers

    path = os.getcwd() + "\\ws"
    mapdl = pyansys.launch_mapdl(run_location=path, interactive_plotting=True, loglevel='ERROR', override=True)
    # """
    # clear existing geometry
    mapdl.finish()
    mapdl.clear()

    th = 50  # [°C]
    tn = 50  # [°C]

    mapdl.prep7()
    # define Material
    for i in range(0, len(materials)):
        mapdl.mptemp(i + 1, th)

        mapdl.mpdata("EX", i + 1, "", materials[i]["e11"])
        mapdl.mpdata("EY", i + 1, "", materials[i]["e22"])
        mapdl.mpdata("EZ", i + 1, "", materials[i]["e33"])

        mapdl.mpdata("PRXY", i + 1, "", materials[i]["pr12"])
        mapdl.mpdata("PRXZ", i + 1, "", materials[i]["pr13"])
        mapdl.mpdata("PRYZ", i + 1, "", materials[i]["pr23"])

        mapdl.mpdata("DENS", i + 1, "", materials[i]["dens"])

        mapdl.mpdata("GXY", i + 1, "", materials[i]["g12"])
        mapdl.mpdata("GXZ", i + 1, "", materials[i]["g13"])
        mapdl.mpdata("GYZ", i + 1, "", materials[i]["g23"])

        # Waermeausdenungsk
        mapdl.mpdata("ALPX", i + 1, "", materials[i]["alp11"])
        mapdl.mpdata("ALPY", i + 1, "", materials[i]["alp22"])
        mapdl.mpdata("ALPZ", i + 1, "", materials[i]["alp33"])

        # Waermeleitfaehigk
        mapdl.mpdata("KXX", i + 1, "", materials[i]["k11"])
        mapdl.mpdata("KYY", i + 1, "", materials[i]["k22"])
        mapdl.mpdata("KZZ", i + 1, "", materials[i]["k33"])

        # Quellkoeffizient
        mapdl.mpdata("BETX", i + 1, "", materials[i]["bet11"])
        mapdl.mpdata("BETY", i + 1, "", materials[i]["bet22"])
        mapdl.mpdata("BETZ", i + 1, "", materials[i]["bet33"])

        # Diffusionskoeffizient
        mapdl.mpdata("DXX", i + 1, "", materials[i]["d11"])
        mapdl.mpdata("DYY", i + 1, "", materials[i]["d22"])
        mapdl.mpdata("DZZ", i + 1, "", materials[i]["d33"])

    mapdl.et(1, "SOLID186")

    mapdl.csys(1)
    mapdl.k(1, str(inner_diameter/2), 0, str(-thickness/2))
    mapdl.k(2, str(inner_diameter/2), 0, str(thickness/2))
    mapdl.k(3, str(outer_diameter/2), 0, str(thickness/2))
    mapdl.k(4, str(outer_diameter/2), 0, str(-thickness/2))
    mapdl.l(1, 2)
    mapdl.l(2, 3)
    mapdl.l(3, 4)
    mapdl.l(4, 1)
    mapdl.al(1, 2, 3, 4)
    mapdl.k(5, 0, 0, str(thickness/2))
    mapdl.k(6, 0, 0, str(-thickness/2))
    mapdl.run("VROT,1,,,,,,5,6")

    arcdiv = 20*4
    darc = 360/arcdiv

    mapdl.lsel("S", "LOC", "X", str(inner_diameter/2))
    mapdl.lsel("R", "LOC", "Z", str(-thickness/2))
    mapdl.lesize("ALL", "", "", str(arcdiv/4))
    mapdl.lsel("S", "LOC", "X", str(inner_diameter/2))
    mapdl.lsel("R", "LOC", "Z", str(thickness/2))
    mapdl.lesize("ALL", "", "", str(arcdiv/4))
    mapdl.lsel("S", "LOC", "X", str(outer_diameter/2))
    mapdl.lsel("R", "LOC", "Z", str(-thickness/2))
    mapdl.lesize("ALL", "", "", str(arcdiv/4))
    mapdl.lsel("S", "LOC", "X", str(outer_diameter/2))
    mapdl.lsel("R", "LOC", "Z", str(thickness/2))
    mapdl.lesize("ALL", "", "", str(arcdiv/4))
    for i in range(0, 4):
        mapdl.lsel("S", "LOC", "Y", str(90*i))
        mapdl.lsel("R", "LOC", "Z", str(-thickness/2))
        mapdl.lesize("ALL", "", "", str(arcdiv/4))
        mapdl.lsel("S", "LOC", "Y", str(90*i))
        mapdl.lsel("R", "LOC", "Z", str(thickness/2))
        mapdl.lesize("ALL", "", "", str(arcdiv/4))
        mapdl.lsel("S", "LOC", "Y", str(90*i))
        mapdl.lesize("ALL", "", "", str(len(angles)))

    mapdl.lsel("ALL")
    mapdl.mshape(0, "3D")
    mapdl.mshkey(1)
    mapdl.vmesh("ALL")
    # mapdl.eplot("ALL")

    for j in range(1, len(angles)+1):
        for i in range(0, 20*4):
            mapdl.csys(1)
            mapdl.esel("S", "CENT", "Y", str(darc*(i + 0.5)))
            mapdl.esel("R", "CENT", "Z", str(thickness*(-1/2+1/len(angles)*(j-1))), str(thickness*(-1/2+1/len(angles)*j)), str(thickness/len(angles)/2))
            mapdl.local(str(j*100+i), 0, 0, 0, 0, str(angles[j-1]+darc*(i + 0.5)))
            mapdl.emodif("all", "ESYS", str(j*100+i))
            if len(angles) == len(materials):
                mapdl.emodif("all", "MAT", str(j))

    mapdl.esel("ALL")

    # Read Stress
    mapdl.run("/SOLU")
    mapdl.csys(1)
    mapdl.omega("", "", omega)
    mapdl.nsel("S", "LOC", "X", inner_diameter / 2)
    mapdl.d("ALL", "UX", 0)
    mapdl.d("ALL", "UY", 0)
    # D,ALL,ROTX,0
    # D,ALL,ROTY,0
    # D,ALL,ROTZ,0
    mapdl.nsel("R", "LOC", "Z", 0)
    mapdl.d("ALL", "UZ", 0)
    mapdl.run("ALLS")
    mapdl.outres("ALL", "ALL")
    mapdl.solve()
    mapdl.finish()

    mapdl.run("/Post1")

    result = mapdl.result
    # result.plot_nodal_solution(0, 'z', show_displacement=True)
    # result.plot_principal_nodal_stress(0, 'S1', show_edges=True, show_axes=True)

    # """
    nodenump, stress = result.principal_nodal_stress(0)
    s1max = np.nanmax(stress[:, -5])
    s2max = np.nanmax(stress[:, -4])

    nodenum, stress = result.nodal_stress(0)
    srmax = np.nanmax(stress[:, -6])
    stmax = np.nanmax(stress[:, -5])

    # Modual-Analysis
    # """
    mapdl.run("/SOLU")
    mapdl.csys(1)
    mapdl.omega("", "", omega)
    mapdl.nsel("S", "LOC", "X", inner_diameter / 2)
    mapdl.d("ALL", "UX", 0)
    mapdl.d("ALL", "UY", 0)
    mapdl.nsel("R", "LOC", "Z", 0)
    mapdl.d("ALL", "UZ", 0)
    mapdl.run("ALLS")
    mapdl.outres("ALL", "ALL")

    mapdl.antype("MODAL")
    mapdl.modopt("lanb", 3)
    mapdl.mxpand("", "", "", "yes")

    mapdl.solve()
    mapdl.finish()
    
    mapdl.run("/Post1")

    result = mapdl.result
    f = result.time_values

    mapdl.exit()
    time.sleep(0.5)

    return s1max, s2max, srmax, stmax, f[0], f[1], f[2]


print(sim_rotor(
    [{
        "e11":      42.5E9  ,
        "e22":      11E9    ,
        "e33":      11E9    ,
        "pr12":     0.28    ,
        "pr13":     0.28    ,
        "pr23":     0.28    ,
        "dens":     1950    ,
        "g12":      4.2E9   ,
        "g13":      4.2E9   ,
        "g23":      2.56E9  ,
        "alp11":    5.7E-6  ,
        "alp22":    45E-6   ,
        "alp33":    45E-6   ,
        "k11":      0.72E3  ,
        "k22":      0.5E3   ,
        "k33":      0.5E3   ,
        "bet11":    0E-3,
        "bet22":    4E-3,
        "bet33":    4E-3,
        "d11":      4.4E-3,
        "d22":      3.1E-3,
        "d33":      3.1E-3
    }], [90, 0, 0, 0, 0, 0]))