Showing direction of surface pressure on element plot

[ay8tb]

Hello,

I'm having some difficulty showing where the surface pressure is on an element plot in order to ensure it's in the right direction and surface. Here is a snippet of the code:

mapdl.geometry.area_select([6,6])
mapdl.sfa('ALL', 1, 'PRES', press)
mapdl.allsel()
mapdl.psf(key=2, kshell=1, color='ON')
mapdl.eplot(vtk=True, cpos="xy", show_edges=True, show_axes=True, plot_bc=True, plot_bc_legend=True)

Any suggestions? I've tried looking at the examples in the documentation but I'm unable to find one that plots before the solution. Thank you!

[germa89]

Hi @ay8tb

It will help a lot if you can share the whole script.

I think probably @mikerife or @pmaroneh can help you better with this. In the meantime, I will use vtk=False to use MAPDL plotting capabilities to get that plot:

mapdl.geometry.area_select([6,6])
mapdl.sfa('ALL', 1, 'PRES', press)
mapdl.allsel()
mapdl.psf(key=2, kshell=1, color='ON')
mapdl.show("png")
mapdl.eplot(vtk=False)
mapdl.show("close")

[ay8tb]

Sure, I've added it below. I use vtk=True since I like to manipulate the GUI and see the boundary conditions, typically I can also see the load placement and direction. I am trying to place the load on the surface that is opposite to the fixed condition in order to simulate an MTS machine test.

# Import necessary libraries
from ansys.mapdl.core import launch_mapdl
from ansys.mapdl.core import launcher

mapdl = launch_mapdl()

# Import ANSYS file
mapdl.clear()
sim_file = mapdl.input(dogbone_geom.anf')
mapdl.finish()
length = 6  # length in inches
width = 0.5  # width in inches
sf_area = 1 * (1 / 8)  # surface area of where force will be applied in inches
sf_load = 1562.50  # load that was applied on MTS
press = sf_load / sf_area

mapdl.allsel()

# properties of coupon
mapdl.prep7()
mapdl.units("BIN")  # U.S. Customary system using inches (in, lbf*s2/in, s, °F).
mapdl.title("Simulation File of Dogbone Coupon")
# define and set materials element type
mapdl.et("1", "SOLID187")  
mapdl.mp("EX", 1, 14100000)  
mapdl.mp("NUXY", 1, 0.311)  
mapdl.mp("DENS", 1, 0.307) 
mapdl.emodif("ALL", "MAT", 1)  

# add meshing
mapdl.smrtsize(1)
mapdl.esize(0.1, 0)
mapdl.mshape(1, "3D")
mapdl.mshkey(0)
mapdl.vmesh("ALL")
# confirm mesh
mapdl.eplot(color="grey", background='b')

mapdl.nsel("S", "LOC", "X", 0, 0)  # fixing in the x direction
mapdl.d("ALL", "UX")
mapdl.d("ALL", "UY")
mapdl.d("ALL", "UZ")

# check boundary conditions initially
mapdl.allsel(mute=True)
# plot the boundary conditions
mapdl.eplot(vtk=True, cpos="xy", show_edges=True, show_axes=True, plot_bc=True, plot_bc_legend=True)

mapdl.aplot(show_area_numbering=True) # used to check what area number to plot surface pressure on

mapdl.geometry.area_select([6,6])
mapdl.sfa('ALL', 1, 'PRES', press)
mapdl.allsel()
mapdl.psf('PRES',key=2, color='ON')
mapdl.show("png")
mapdl.eplot(vtk=False)
mapdl.show("close")

# run solution
mapdl.run("/SOL")
mapdl.antype(0)
mapdl.outres("ALL", "ALL")
mapdl.run("/STATUS,SOLU")
sol_output = mapdl.solve()

[germa89]

We are going to need the dogbone_geom.anf file.

[ay8tb]

Sorry about that @germa89, I've attached it now.
dogbone_geom.zip

[mikerife]

Hi @ay8tb
The PyMAPDL plotting does not yet support pressure symbols - please see the following for what all bc/load symbols it does currently support:

https://mapdl.docs.pyansys.com/version/stable/api/_autosummary/ansys.mapdl.core.plotting.general_plotter.html

In the PyMAPDL script you provided since the pressure is applied to an area the pressure symbols do not show up when plotting elements via the:

mapdl.sfa('ALL', 1, 'PRES', press)
mapdl.allsel()
mapdl.psf('PRES',key=2, color='ON')
mapdl.show("png")
mapdl.eplot(vtk=False)
mapdl.show("close")

You may want to either replace the eplot with a aplot(vtk=False) or use the open_gui() to open the model in MAPDL interactive and double check that the model is set up correctly. Or you can issue a SFALIST command to list the applied surface pressures.

Mike

[ay8tb]

Hi @mikerife,

I'm getting back to this after a while and I've tried to do a few things differently:

1. Instead of pressure, use the node method found here: https://mapdl.docs.pyansys.com/version/stable/examples/gallery_examples/00-mapdl-examples/3d_notch.html
2. Instead of loading in pyMAPDL, using the open_gui() method to place the load on the nodes and then continue analysis in pyMAPDL.
3. Use the log_apdl function to record input from (1) and compare with Workbench. This method didn't work at all (geometry wasn't imported, loads weren't imported, was confusing to navigate.

Both of these methods give me similar results to the Workbench solution (and the one that is equivalent to hand calculations) but the results are still a bit off. I've checked that the material properties and load magnitude is equivalent in both methods. Is there anything in the script that doesn't match a fixed boundary condition at x=0 and a positive load at x=length of 4600 N? In the Workbench code, I applied the 4.6 kN load to the face of the structure, would that make a difference?

#Import necessary libraries
from ansys.mapdl.core import launch_mapdl
from ansys.mapdl import core as pymapdl
from ansys.mapdl.core import launcher
import numpy as np
import pyvista as pv
import os
import pandas as pd

print(launcher.get_default_ansys_path())  # making sure anysy path is the same as where it's stored
pymapdl.change_default_ansys_path('C:\\Program Files\\ANSYS2\\ANSYS Inc\\v222\\ansys\\bin\\winx64\\ANSYS222.exe')
mapdl = launch_mapdl(start_timeout=300)
result_path = "/"  # change on where Python files are, esp the .rst file

# Import ANSYS file 
mapdl.clear()
sim_file = mapdl.input('C:\\Users\\XXXX\\XXX\\Desktop\\Geom_trial_mm.anf')
mapdl.finish()
mm_converter = 25.4  # to convert base in inches to mm to match DIC output
#length = 6 * mm_converter  # length in inches
length = 152.4 # lenght in mm
width = 1 * mm_converter  # width in inches
half_width = width/2

mapdl.allsel()

# properties of coupon
mapdl.prep7()
mapdl.units("MPA")  # International system using mm (mm, Mg, s, °C)
mapdl.title("Simulation File of Dogbone Coupon")

elem_type = 'SOLID185'
mapdl.et(1, "SOLID185")  # SOLID185 - 3-D 8-Node Hexahedral Structural Solid
mapdl.mp("EX", 1, 211400)  # elastic modulus, element x direction in MPa
mapdl.mp("NUXY", 1, 0.29)  # Major Poisson's ratio, x-y plane - CAN CHANGE THIS
mapdl.mp("DENS", 1, 0.000000000785)  # mass density in tonne/mm3
mapdl.emodif("ALL", "MAT", 1)  # defining material and modifying

#
mapdl.esize(0.4)
mapdl.et(1,185)
mapdl.vsweep(  "ALL")  # Fills an existing unmeshed volume with elements by sweeping mesh from an adjacent area through volume.



print("Completed mesh")
# fix both ends of coupon
# first, fix at the 0 location
# Fix the left-hand side.

plotter = pv.Plotter(shape=(1, 3))

plotter.subplot(0, 0)
mapdl.lsla("S")
mapdl.lplot(vtk=True, show_keypoint_numbering=True, plotter=plotter)
mapdl.lsel("all")
mapdl.nsel("S", "LOC", "X", 0)
selected_nodes = mapdl.mesh.nnum
print(selected_nodes)
mapdl.d("ALL", "UX")

# Fix a few nodes on the left-hand side of the plate in the Y and Z
# direction.  Otherwise, the mesh would be allowed to move in the y
# direction and would be an improperly constrained mesh.
mapdl.nsel("R", "LOC", "Y", half_width)
selected_nodes = mapdl.mesh.nnum
print(selected_nodes)
mapdl.d("ALL", "UY")
mapdl.d("ALL", "UZ")

#mapdl.nsel("S", "LOC", "X", length)

# Check if any nodes are selected at the specified location
#selected_nodes = mapdl.mesh.nnum
#print(selected_nodes)

# Verify that only the nodes at length have been selected:
# assert np.unique(mapdl.mesh.nodes[:, 0]) == length

# Next, couple the DOF for these nodes.  This lets us provide a force
# to one node that will be spread throughout all nodes in this coupled
# set.
mapdl.cp(5, "UX", "ALL")

# Select a single node in this set and apply a force to it
# We use "R" to re-select from the current node group
mapdl.nsel("R", "LOC", "Y", half_width)  # selects more than one
selected_nodes = mapdl.mesh.nnum
single_node = mapdl.mesh.nnum[0]
mapdl.nsel("S", "NODE", vmin=single_node, vmax=single_node)
mapdl.f("ALL", "FX", 4600)

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


# run solution
mapdl.run("/SOL")
mapdl.antype(0)
mapdl.outres("ALL", "ALL")
mapdl.run("/STATUS,SOLU")
sol_output = mapdl.solve()

# grab the result from the ``mapdl`` instance
result = mapdl.result