mapdl.non_interactive crashes randomly

[germa89]

=== Reported internally DK ===

It seems that running many times with the mapdl.non_interactive generates a MapdlExitedError error.

CRITICAL -  -  logging - handle_exception - Uncaught exception
Traceback (most recent call last):
  File "C:/Work/E_Submodel/Issue_LoadSets/2d_plate_with_a_hole_new.py", line 184, in <module>
    stress_con = compute_stress_con(ratio)
  File "C:/Work/E_Submodel/Issue_LoadSets/2d_plate_with_a_hole_new.py", line 72, in compute_stress_con
    mapdl.emodif("ALL", "MAT", 1)
  File "C:\ENERCON\Python38\lib\site-packages\ansys\mapdl\core\mapdl.py", line 512, in __exit__
    self._parent()._flush_stored()
  File "C:\ENERCON\Python38\lib\site-packages\ansys\mapdl\core\mapdl_grpc.py", line 1412, in _flush_stored
    out = self.input(
  File "C:\ENERCON\Python38\lib\site-packages\ansys\mapdl\core\errors.py", line 131, in wrapper
    raise MapdlExitedError("MAPDL server connection terminated") from None
ansys.mapdl.core.errors.MapdlExitedError: MAPDL server connection terminated

Process finished with exit code 1

Code:

"""
.. _ref_plane_stress_concentration:

MAPDL 2D Plane Stress Concentration Analysis
--------------------------------------------

"""

import matplotlib.pyplot as plt
import numpy as np

from ansys.mapdl.core import launch_mapdl

mapdl = launch_mapdl(run_location=r'C:\Work\E_Submodel\20220509\temp', additional_switches= '-smp', override=True,
                     loglevel="DEBUG",  # to log everything
                     log_apdl=r"C:\Work\E_Submodel\20220509\my_apdl_log.log"  # to log the MAPDL commands
                     )

length = 0.4
width = 0.1

ratio = 0.3  # diameter/width
diameter = width * ratio
radius = diameter * 0.5

def compute_stress_con(ratio):
    """Compute the stress concentration for plate with a hole loaded
    with a uniaxial force.
    """
    with mapdl.non_interactive:

        mapdl.finish()
        mapdl.clear("NOSTART")
        mapdl.prep7()
        mapdl.units("SI")  # SI - International system (m, kg, s, K).

        # define a PLANE183 element type with thickness
        mapdl.et(1, "PLANE183", kop3=3)
        mapdl.r(1, 0.001)  # thickness of 0.001 meters)

        # Define a material (nominal steel in SI)
        mapdl.mp("EX", 1, 210e9)  # Elastic moduli in Pa (kg/(m*s**2))
        mapdl.mp("DENS", 1, 7800)  # Density in kg/m3
        mapdl.mp("NUXY", 1, 0.3)  # Poisson's Ratio
        mapdl.emodif("ALL", "MAT", 1)

    # Geometry
    # ~~~~~~~~
    # Create a rectangular area with the hole in the middle
    diameter = width * ratio
    radius = diameter * 0.5

    # create the rectangle
    rect_anum = mapdl.blc4(width=length, height=width)

    # create a circle in the middle of the rectangle
    circ_anum = mapdl.cyl4(length / 2, width / 2, radius)

    # Note how pyansys parses the output and returns the area numbers
    # created by each command.  This can be used to execute a boolean
    # operation on these areas to cut the circle out of the rectangle.
    plate_with_hole_anum = mapdl.asba(rect_anum, circ_anum)

    with mapdl.non_interactive:
        # Meshing
        # ~~~~~~~
        # Mesh the plate using a higher density near the hole and a lower
        # density for the remainder of the plate

        mapdl.aclear("all")

        # ensure there are at least 100 elements around the hole
        hole_esize = np.pi * diameter / 100  # 0.0002
        plate_esize = 0.001

        # increased the density of the mesh at the center
        mapdl.lsel("S", "LINE", vmin=5, vmax=8)
        mapdl.lesize("ALL", hole_esize, kforc=1)
        mapdl.lsel("ALL")

        # Decrease the area mesh expansion.  This ensures that the mesh
        # remains fine nearby the hole
        mapdl.mopt("EXPND", 0.7)  # default 1

        mapdl.esize(plate_esize)
        mapdl.amesh(plate_with_hole_anum)

        ###############################################################################
        # Boundary Conditions
        # ~~~~~~~~~~~~~~~~~~~
        # Fix the left-hand side of the plate in the X direction
        mapdl.nsel("S", "LOC", "X", 0)
        mapdl.d("ALL", "UX")

        # Fix a single node on the left-hand side of the plate in the Y direction
        mapdl.nsel("R", "LOC", "Y", width / 2)

    assert mapdl.mesh.n_node == 1

    with mapdl.non_interactive:
        mapdl.d("ALL", "UY")

        # Apply a force on the right-hand side of the plate.  For this
        # example, we select the right-hand side of the plate.
        mapdl.nsel("S", "LOC", "X", length)

        # Next, couple the DOF for these nodes
        mapdl.cp(5, "UX", "ALL")

        # Again, select a single node in this set and apply a force to it
        mapdl.nsel("r", "loc", "y", width / 2)
        mapdl.f("ALL", "FX", 1000)

        # finally, be sure to select all nodes again to solve the entire solution
        mapdl.allsel()

        # Solve the Static Problem
        # ~~~~~~~~~~~~~~~~~~~~~~~~
        mapdl.run("/SOLU")
        mapdl.antype("STATIC")
        mapdl.solve()

    return 1
    # Post-Processing
    # ~~~~~~~~~~~~~~~
    # grab the stress from the result
    result = mapdl.result
    nnum, stress = result.principal_nodal_stress(0)
    von_mises = stress[:, -1]
    max_stress = np.nanmax(von_mises)

    # compare to the "far field" stress by getting the mean value of the
    # stress at the wall
    mask = result.mesh.nodes[:, 0] == length
    far_field_stress = np.nanmean(von_mises[mask])

    # adjust by the cross sectional area at the hole
    adj = width / (width - diameter)
    stress_adj = far_field_stress * adj

    # finally, compute the stress concentration
    return max_stress / stress_adj


k_t_exp = []
ratios = np.linspace(0.01, 0.5, 60)
print("    Ratio  : Stress Concentration (K_t)")
import time

t1 = time.time()
num_of_jobs = len(ratios)
for i, ratio in enumerate(ratios, start=1):
    print(f"job {i} of {num_of_jobs}")
    t2 = time.time()
    stress_con = compute_stress_con(ratio)
    print("%10.4f : %10.4f" % (ratio, stress_con))
    k_t_exp.append(stress_con)
    t3 = time.time()
    print("time one job: ", t3-t2)
    print("total time: ", t3-t1)


mapdl.exit()

[lschuler1]

Hi @germa89,
I think I observed the same behavior (mapdl crashing randomly) but without using mapdl.non_interactive.
I am using mapdl version 0.61.3, and for instance I get this error when running a code that was working:

[germa89]

mmhh... without using non_interactive?? ... If you can, can you upload your code here, so I can have a look??

[lschuler1]

Hi @germa89 ,

I cannot give you my entire code, but here is the part that caused the issue. It's a loop to obtain the right hand side vector of a thermal problem for different heat generation values. MAPDL crashed randomly at different iterations of the loop when calling mapdl.math.rhs().

for i in range(len(pJ[0])): # for each mode of the gradient prod.
    # Print mode number
    print('Compute thermal force number %i of %i'%(i + 1,len(pJ[0])))

    # Compute the thermal force
    mapdl.slashsolu()
    mapdl.wrfull('1')
    mapdl.kuse('1')
    mapdl.time(1) ; mapdl.autots('OFF') ; mapdl.nsubst(1)
    mapdl.outres('ERASE')
    mapdl.outres('ALL','NONE') # no data written to data files

    # Clear previous load
    mapdl.bfedele('ALL','HGEN')

    # Sets uniform initial temperature to 0 and deletes convection
    # otherwise the RHS is impacted
    mapdl.tunif('0')
    mapdl.sfdele('ALL','CONV')

    # Create commands to apply the load
    cmd = ['bfeblock,3,hgen,{Emax},{nE},0'
           .format(Emax=len(elem), nE=len(elem))]
    # 3: nb. of fiels, nE: maximum element number defined, number of selected elements, 0: non tabular input
    cmd.append('(2i9,e20.9e3)')
    for e in range(len(elem)):
        cmd.append('%9i%9i%20.9e' %(e+1,1,pJ[0][i][e]))

    cmd.append('bfe,end,loc,-1,')

    # Pass it to MAPDL and solve
    mapdl.input_strings(cmd)
    mapdl.finish()

    # Solve
    mapdl.slashsolu()
    mapdl.solve()
    mapdl.finish()
    
    # combine .full files from the processors
    mapdl.aux2()
    mapdl.combine('full')
    mapdl.finish()

    F = mapdl.math.rhs(asarray=True)
    Fth[0][i] = F.reshape([len(F),1])

[germa89]

I will have a look in the coming months, I'm a bit caught with other stuff. Please feel free to ping me if no progress is made in the next 2 weeks.

[lschuler1]

No problem. I'haven't encountered the issue for a while so there is no hurry.