How to extract Element Nodal Forces with PyDPF

Options

Hello,

I am currently trying to develop a script that uses the elementary nodal forces and I have to perform the extraction with PyDPF from the result file. I am working with Ansys 2023R1 and here is my current script:

import os
from ansys.dpf import core as dpf

input_file = "file.rst"
bolt_selection = "SELECTION"

axis = 0
coord_sys = 0

# Select the ansys version 2023 for DPF
ansys_path = os.getenv("AWP_ROOT232")
if ansys_path is None:
    ansys_path = os.getenv("AWP_ROOT231")
# print(ansys_path)

# Start the DPF server
server = dpf.start_local_server(ansys_path=ansys_path)
# print(server)

# Open the model
model = dpf.Model(input_file)
# print(model)

# Read the results
results = model.results
print(results)

# Extract the scoping of the named selection
mesh_scoping = model.metadata.named_selection(bolt_selection)
print(mesh_scoping)
print(mesh_scoping.ids)

# Extract nodal forces on selection
enf = results.element_nodal_forces(
    time_scoping=1, mesh_scoping=mesh_scoping
)
fields = enf.eval()
print(fields)
field = fields[0]
print(field.data)

The script runs, I will attach to this request an archive containing the model to provide the data, but my issue is that I do not obtain the same results with ENF in Ansys. In this specific test case, I am expecting to obtain -25N in x direction but the script returns zero.

After investigating I think this is related to the scoping I am using. Element nodal forces seems to be defined on Element Nodal location while my scoping coming from the named selection is defined on Nodal location. However I don't know how to fix this to get the same results as the one obtained with Mechanical.

Would you know how to proceed?

Thank you in advance!

Tagged:

Answers

  • Paul Creusy
    Paul Creusy Member Posts: 5
    Name Dropper First Comment
    Options

    I have tried to develop a script to realise the average on the nodes that are concerned by the selection for each adjacent element and it works perfectly for simple cases but I am still unable to recover the same results as in Ansys.

    def get_nodes_hit(nodes):
        nodes_hit = {}
        for node in nodes:
            nodes_hit[node] = {
                "elements": [],
                "local_ids": [],
            }
        for el in model.metadata.meshed_region.elements:
            el_nodes = el.node_ids
            for i, node in enumerate(nodes):
                for j in range(len(el_nodes)):
                    if node == el_nodes[j]:
                        if el.type in accepted_types:
                            nodes_hit[node]["elements"].append(el.id)
                            nodes_hit[node]["local_ids"].append(j + 1)
        return nodes_hit
    
    
    def get_results(nodes, time_scoping, axis):
        # Get the nodes hit to know on which elements to extract
        nodes_hit = get_nodes_hit(nodes)
    
        # Prepare the list of elements for which an extraction is required
        scoping_elements = []
        for node in nodes:
            for el in nodes_hit[node]["elements"]:
                if el not in scoping_elements:
                    scoping_elements.append(el)
    
        # Create the scoping
        scoping = mesh_scoping_factory.elemental_scoping(scoping_elements)
    
        # Get the results on the elements
        el_res = results.element_nodal_forces(
            time_scoping=time_scoping, mesh_scoping=scoping
        )
        fields = el_res.eval()
        field = fields[0]
    
        # Convert element results to nodal results
        node_res = []
        for node in nodes:
            temp = 0
            for i in range(len(nodes_hit[node]["elements"])):
                current_element = nodes_hit[node]["elements"][i]
                current_local_id = nodes_hit[node]["local_ids"][i]
                for j, el in enumerate(scoping_elements):
                    if current_element == el:
                        temp+= field.get_entity_data_by_id(el)[current_local_id - 1, axis]
            node_res.append(temp)
    
        return node_res
    
  • Paul Creusy
    Paul Creusy Member Posts: 5
    Name Dropper First Comment
    Options

    I have tried to develop a new script to compute manually the average on the nodes and it works perfectly for simple results but it fails for more complex models. Would you know how to solve it?

    def get_nodes_hit(nodes):
        nodes_hit = {}
        for node in nodes:
            nodes_hit[node] = {
                "elements": [],
                "local_ids": [],
            }
        for el in model.metadata.meshed_region.elements:
            el_nodes = el.node_ids
            for i, node in enumerate(nodes):
                for j in range(len(el_nodes)):
                    if node == el_nodes[j]:
                        if el.type in accepted_types:
                            nodes_hit[node]["elements"].append(el.id)
                            nodes_hit[node]["local_ids"].append(j + 1)
        return nodes_hit
    
    
    def get_results(nodes, time_scoping, axis):
        # Get the nodes hit to know on which elements to extract
        nodes_hit = get_nodes_hit(nodes)
    
        # Prepare the list of elements for which an extraction is required
        scoping_elements = []
        for node in nodes:
            for el in nodes_hit[node]["elements"]:
                if el not in scoping_elements:
                    scoping_elements.append(el)
    
        # Create the scoping
        scoping = mesh_scoping_factory.elemental_scoping(scoping_elements)
    
        # Get the results on the elements
        el_res = results.element_nodal_forces(
            time_scoping=time_scoping, mesh_scoping=scoping
        )
        fields = el_res.eval()
        field = fields[0]
    
        # Convert element results to nodal results
        node_res = []
        for node in nodes:
            temp = 0
            for i in range(len(nodes_hit[node]["elements"])):
                current_element = nodes_hit[node]["elements"][i]
                current_local_id = nodes_hit[node]["local_ids"][i]
                for j, el in enumerate(scoping_elements):
                    if current_element == el:
                        temp += field.get_entity_data_by_id(el)[current_local_id - 1, axis]
            node_res.append(temp)
    
        return node_res
    
    
  • Paul Creusy
    Paul Creusy Member Posts: 5
    Name Dropper First Comment
    Options

    I also tried another method using DPF directly within Mechanical but then again I am unable to obtain the same results as the one displayed with ENFO.

    Here is the script that I used:

    import mech_dpf
    
    import Ans.DataProcessing as dpf
    
    import os
    
    
    
    mech_dpf.setExtAPI(ExtAPI)
    
    analyses = ExtAPI.DataModel.AnalysisList
    
    analysis = analyses[0]
    
    analysis_dir_path = analysis.WorkingDir
    
    path = os.path.join(analysis_dir_path, 'file.rst')
    
    model = dpf.Model(path)
    
    # Define mesh scoping
    nsName_nodal = model.AvailableNamedSelections[1]
    print(nsName_nodal)
    mesh_scoping = dpf.MeshScopingFactory.NamedSelectionScoping(nsName_nodal, model)
    mesh_scoping.Location = 'Nodal'
    
    # Define time scoping
    timeScop = dpf.Scoping()
    timeScop.Ids = [289]
    
    
    dataSources = mech_dpf.GetDataSources(0)
    #res = dpf.operators.result.stress_X()
    res = dpf.operators.result.element_nodal_forces()
    
    res.inputs.data_sources.Connect(dataSources)
    
    res.inputs.mesh_scoping.Connect(mesh_scoping)
    res.inputs.time_scoping.Connect(timeScop)
    
    fields = res.outputs.fields_container.GetData()
    
    print(fields)
    real_part = fields[0].Data
    
    real_part = real_part.ToArray()
    print(real_part)
    
    enfoz = []
    for i in range(2,390,3):
        enfoz.append(real_part[i])
    
    print(min(enfoz))
    print(max(enfoz))
    print(sum(enfoz))
    

    I don't see many other options anymore, I would really appreciate if you have some ideas on how to solve it!
    Thank you in advance!

  • Pierre Thieffry
    Pierre Thieffry Member, Moderator, Employee Posts: 101
    First Anniversary Ansys Employee Solution Developer Community of Practice Member Photogenic
    Options

    @Paul Creusy did you try dpf.operators.result.raw_reaction_force ? This may lead to better results.

  • Paul Creusy
    Paul Creusy Member Posts: 5
    Name Dropper First Comment
    Options

    @Pierre Thieffry thank you for your answer.
    I guess that you mean dpf.operators.result.reaction_force because I obtain an error message if I add the raw. When I try to access those results, I obtain only zeros on my surface, even if I enabled the output control for reaction forces.