How to extract Element Nodal Forces with PyDPF

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: 6
    Name Dropper First Comment
    **

    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: 6
    Name Dropper First Comment
    **

    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: 6
    Name Dropper First Comment
    **

    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: 107
    25 Answers 10 Comments 25 Likes First Anniversary
    ✭✭✭✭

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

  • Paul Creusy
    Paul Creusy Member Posts: 6
    Name Dropper First Comment
    **

    @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.

  • Pierre Thieffry
    Pierre Thieffry Member, Moderator, Employee Posts: 107
    25 Answers 10 Comments 25 Likes First Anniversary
    ✭✭✭✭
    edited July 11

    @Paul Creusy Can you try this? I see consistent results with ENFO. I only scoped to solid elements, then rescoped to your selection

    import os
    from ansys.dpf import core as dpf
    
    input_file = r"D:\Support\forum\nodal forces\2_el_beam_files\dp0\SYS\MECH\file.rst"
    bolt_selection = "SELECTION"
    
    axis = 0
    coord_sys = 0
    
    
    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)
    # Start the DPF server
    #server = dpf.start_local_server()
    print(server)
    
    ds=dpf.DataSources(input_file)
    # Open the model
    model = dpf.Model(ds)
    mesh=model.metadata.meshed_region
    # 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)
    
    
    elem_scoping=dpf.operators.scoping.on_mesh_property(requested_location="Elemental",
                                                        property_name='apdl_element_type',
                                                        property_id=185,mesh=mesh)
    
    # Extract nodal forces on selection
    enf = dpf.operators.result.element_nodal_forces(data_sources=ds,
                                                    time_scoping=1,                                                 
                                                    requested_location='Nodal',
                                                    mesh_scoping=elem_scoping)
    
    scoped_enf=dpf.operators.scoping.rescope(enf.outputs.fields_container()[0],
                                             mesh_scoping=mesh_scoping)
    
    all_forces=enf.eval()
    fields = scoped_enf.outputs.fields_as_field()
    print(fields)
    print(fields.data)
    
  • Paul Creusy
    Paul Creusy Member Posts: 6
    Name Dropper First Comment
    **

    Hello @Pierre Thieffry ,
    Thank you very much for your answer, the extraction script you propose provides the expected results for me as well.
    However I want to realize the extraction over all the frequencies available for a model and the rescoping seems to contain only one time_scoping. Therefore I have included in my script the lines 43 and 48 in a loop to select a new time scoping at each iteration, but this method is very slow, each iteration taking around 12 secs, most of the time being spent on the rescoping.

    Would you know a way to realize the rescoping only once and then extract the values for the different frequencies one by one ?

    Thank you again for your help!