Creating .frd result files using pypdf

M · December 2023

A recent request from an Ansys user to create results in a .frd format.
This is possible using the pyansys pydpf module. PyDPF
There are some minor issues with pyramid elements, Ansys will generate these in some transition zones and they are not used in the .frd documentation I found. The solution is to convert pyramids into tets. One very minor drawback is that quadratic pyramids will lose the mid side node results. Anyone who comes up with a smart solution is welcome to add to this code.
Note the wedges have not been tested but there is a line early in the code (etypedict) to set the conversion. But wedges should work.
Please read the Ansys documentation for SOLID185/186 for more information about the nodes, connectivity and element information:
solid185
solid186

Note the script was tested on beams, lines, shells and solids but not thoroughly. Check your results carefully.
The assumption is a .rst (structural) file from an Ansys analysis is available but dpf can also read .rth (thermal) files so it can be adapted. Currently the script requires displacement results but will try to find stress results.
The code is working and can be improved for both timing (a timer is there if you want to test) and readability.
The .frd translation was possible using: CalculiX

The ASCII file that is generated requires a certain syntax for nodes/elements/results that is described in the document above.

For more on pydpf please visit the pyansys pydpf page:
PyDPF

pydpf Ansys result

same result viewed from .frd file

M · December 2023

"""
Demo script to read in a .rst file and create a .frd ASCII result file (calculix readable)
"""
timer = True

from ansys.dpf import core as dpf
import os, webbrowser, easygui, time

# file data
export_file_type = 'frd'
default_search_directory = 'D:\\'
result_file_type = 'rst'
here = os.getcwd()

# Element set up for conversion from Ansys to Cclx
pyramids = False
etypedict ={'Line2':11,'Quad4':9,'Quad8':10,'Hex8':1,'Hex20':4, 'Tri6':2, 'Tet4':3,'Tet10':6, 'Pyramid5': 3, 'Pyramid13':3} 
# yet untested: prism options, hard to generate: prism6 and prism15? In theory they should work 1 to 1 w/ calculix
etlist = list(etypedict.values())
et_exceptions = [3,4] # element types that need special handling
# ###############################################################################
# Create a model targeting a file.rst result file
path2file = easygui.fileopenbox(msg="Choose a file", default=default_search_directory+'*.' +result_file_type ,filetypes='*.' + result_file_type)

file_name = path2file.split(os.sep)[-1].split('.')[0]
export_file = os.path.join(here,file_name + '.' + export_file_type)
model = dpf.Model(path2file)

metadata = model.metadata
antype = metadata.result_info.analysis_type
if antype == 'static':
    antype = 1
else:
    antype = 0
n_steps = metadata.time_freq_support.n_sets
mesh = model.metadata.meshed_region

# model data
nodes = model.metadata.meshed_region.nodes.scoping.ids 
elements = mesh.elements.scoping.ids
elements = sorted(elements)
materials = mesh.property_field('mat').data
element_scoping = model.metadata.meshed_region.elements.scoping

def gen_displacement_string(data,nodes,nset,antype):
    fmt1 =  '{:12.5E}'
    nnum = len(nodes)
    #t = nset+1
    col = len(data.data[0])
    preamble = '{:2}{:3d}C{:6d}{:12.5E}{:12d}{:20}{:2d}{:5d}{:10}{:2d}\n'.format(' ',100,100+nset,nset,nnum,' ',antype, nset, ' ',1) 
    preamble = preamble +'{:1}{:2d}{:2}{:6s}{:7d}{:5d}\n'.format(' ', -4, ' ', 'DISP',col+1,1)
    ctr = 1
    for c in range(1,col+1):
        preamble = preamble +'{:1}{:2d}{:2}{:8s}{:5d}{:5d}{:5d}{:5d}\n'.format(' ',-5,' ', 'D'+str(c),1,2,ctr,0)
        ctr +=1
    preamble = preamble + '{:1}{:2d}{:2}{:8s}{:5d}{:5d}{:5d}{:5d}{:5d}{:8s}     \n'.format(' ',-5,' ','ALL',1,2,0,0,1,'ALL')
    node_u = []
    for n in nodes:
        node_u.append([n,*data.get_entity_data_by_id(n)[0]])
    node_u = [i for j in node_u for i in j]
    fmt = (' -1{:10.0f}'+col*fmt1 + '\n')*nnum
    disp = fmt.format(*node_u)
    disp = preamble + disp + ' -3\n'
    return disp


def gen_stress_string(stress_op,nodes,nset,antype):
    """
    Get the stress results, return string for writing to file
    """
    stress_op.inputs.time_scoping.connect([nset])
    stress_fields = stress_op.outputs.fields_container()
    stress_fields_nodal_op=dpf.operators.averaging.elemental_nodal_to_nodal(field=stress_fields[0])
    stress = stress_fields_nodal_op.eval()
    node_s = []

    if stress.size == 0:
        return 0
    else:
        ss =['SXX','SYY','SZZ','SXY','SYZ','SZX']
        count1 = [1,2,3,1,2,3]
        count2 = [1,2,3,2,3,1]
        nnum = len(nodes)
        col = len(stress.data[0])
        preamble = '{:2}{:2d}C{:6d}{:12.5E}{:12d}{:20}{:2d}{:5d}{:10}{:2d}\n'.format(' ',100,100+nset,nset,nnum,' ',antype, nset, ' ',1) 
        preamble = preamble +' {:2d}  {:8s}{:5d}{:5d}\n'.format(-4,'STRESS',col,1)

        for s,c1,c2 in zip(ss,count1,count2):
            preamble = preamble + ' {:2d}  {:5}{:8d}{:5d}{:5d}{:5d}\n'.format(-5,s,1,4,c1,c2)

        for n in nodes:
            node_s.append([n,*stress.get_entity_data_by_id(n)[0]])
        node_s = [i for j in node_s for i in j]
        fmt = (' -1{:10.0f}' + col*'{:12.5E}' + '\n')*nnum
        node_s = fmt.format(*node_s)
        node_s = preamble + node_s + ' -3\n'
        return node_s

# # # Get the results
displacements = model.results.displacement
stress_operator = dpf.operators.result.stress()
stress_operator.inputs.data_sources.connect(model)
# stress test y/n
try:
    stress_operator.inputs.time_scoping.connect([1])
    stress_fields = stress_operator.outputs.fields_container()
    stress_fields_nodal_op=dpf.operators.averaging.elemental_nodal_to_nodal(field=stress_fields[0])
    stress = stress_fields_nodal_op.eval()
    if stress.size == 0:
        stress_avail = False
    else:
        stress_avail = True
except:
    stress_avail = False
    pass

t1 = time.time()
elist = []
fmt = []

fmt1 = '{:10.0f}'
fmt2 = '{:5.0f}'
# Generate formatted element block as string
element_nodes = []
for element,material in zip(elements,materials):
    edata = mesh.elements.element_by_id(element)
    node_indx = edata.connectivity
    try:
        etype = etypedict[edata.type.name]
    except:
        etype = 1
    elist.append([element,etype,0,material])
    element_conx = nodes[node_indx]
    #print(element, edata.type.name,etype,len(element_conx))
    n_e_nodes = len(element_conx)
    if n_e_nodes > 1 : # write element number and type to file
        fmt.append(' -1' + fmt1 + 3*fmt2 + '\n')
        element_nodes.extend(element_conx)
        if etype not in et_exceptions:#  [1,]:#n_e_nodes <= 10:# and etype != 33:
            elist.append(element_conx)
            fmt.append(' -2'+fmt1*n_e_nodes + '\n')
        elif etype == 4:
            elist.append(element_conx[0:10])
            fmt.append(' -2'+fmt1*10 + '\n')
            node_swap_list =[10, 11, 16, 17, 18, 19, 12, 13, 14, 15]
            ordered_list = [element_conx[i] for i in node_swap_list]
            elist.append(ordered_list) 
            fmt.append(' -2'+fmt1*10 + '\n')

        elif etype == 3:# or etype == 6: 
            if len(element_conx) == 4: # true tets 
                elist.append(element_conx)
                fmt.append(' -2'+fmt1*4 + '\n')
            if len(element_conx) == 5 or len(element_conx) == 13:# Ansys 13 node pyramids have to be converted to 2 type 4 (mid side nodes ignored)
                if not pyramids:
                    print('4 node pyramids being added to element list')
                    pyr5_2_elist = {}
                if len(element_conx) == 13:
                    print('13 node pyramids being converted to 4 node pyramids and added to element list')
                pyramids = True
                elist.append([element_conx[i] for i in [0,1,2,4]])
                fmt.append(' -2'+fmt1*4 + '\n')
                pyr5_2_elist[element] = {'mat': material, 'etype':etype, 'nodes' : [element_conx[i] for i in [2,3,0,4]], 'formatting':(' -2'+fmt1*4 + '\n')}

ctr = len(elements)+1
# here any extra generated pyramid elements can be added to the element list
if pyramids:
    for pyr_elem in pyr5_2_elist:
        pyramid_element_data = pyr5_2_elist[pyr_elem]
        node_l = pyramid_element_data['nodes']
        elist.append([ctr,pyramid_element_data['etype'],0,pyramid_element_data['mat']])
        fmt.append(' -1' + fmt1 + 3*fmt2 + '\n')
        elist.append(node_l)
        fmt.append(' -2'+fmt1*len(node_l) + '\n')
        ctr+=1

ecount = ctr-1    
elem_data = [i for j in elist for i in j]

fmt_data = ''.join([i for i in fmt])
elem_data = fmt_data.format(*elem_data)

# element data to write
element_write_data = '    {:2s}{:18s}{:12d}{:37s}{:1d}\n'.format('3C',' ',ecount,' ',1)
element_write_data = element_write_data + elem_data + ' {:2d}\n'.format(-3)

# Node data
t2 = time.time()
nodes = set(element_nodes)
dt1 = str(round(t2-t1,2))
num_elem = str(len(elements))
if timer: print(f'Time to format {num_elem} elements data : {dt1} seconds.') 

fmtNode = '{: 12.5e}'
n_data = []
for n in nodes:
    loc = mesh.nodes.node_by_id(n).coordinates
    n_data.append([n,*loc])
fmt = (' -1{:10.0f}'+3*fmtNode + '\n')*len(n_data)
node_data = [i for j in n_data for i in j]
node_data = fmt.format(*node_data) 
nnum = len(nodes)
node_write_data = '    {:2s}{:18s}{:12d}{:37s}{:1d}\n'.format('2C',' ',nnum ,' ',1)
node_write_data = node_write_data + node_data + ' {:2d}\n'.format(-3)

M · December 2023

t3 = time.time()
dt1 = str(round(t3-t2,2))
num_elem = str(len(elements))
if timer: print(f'time to get node data: {dt1} seconds.') 
with open(export_file,'w',newline='') as f:
    f.write('    1C\n')
    f.write('    1UModel:' + str(file_name) + '\n')
    f.write(node_write_data)
    f.write(element_write_data)

    ctr = 1
    for step in range(1,n_steps+1):
        tx =time.time()
        disp = displacements.on_time_scoping(step).eval()[0]
        disp_block = gen_displacement_string(disp,nodes,step,antype)
        f.write(disp_block)
        tx1 = time.time()
        if timer: print('Time to write displacement block set {:3.0f}: {:3.1f} seconds'.format(step,round(tx1 - tx,2)))
        # stress block
        #line = '    1PSTEP%26s%12s%12s\n' % (ctr,1,step) # for writing parameters to file?
        #f.write(line)
        if stress_avail:
            stress_block = gen_stress_string(stress_operator,nodes,step,antype)
            if stress_block != 0:
                f.writelines(stress_block)
            tx2 = time.time()
            if timer: print('Time to write stress block set {:3.0f}: {:3.1f} seconds'.format(step,round(tx2 - tx1,2)))
    f.write(' 9999\n')        

t4 = time.time()
print('Time to write data for %10s elements to file: %5s seconds' % (len(elements), str(round(t4-t3,2))))  
print('Overall time : %5.6s seconds' % round(t4-t1,2))  

#webbrowser.open(export_file) # fancify, opens up the .frd with your default frd opener
model.metadata.release_streams() # release the result file (.rst)

M · December 2023

Apologies for having the break the script into two parts, it ended up being longer than acceptable to the web interface here.

Creating .frd result files using pypdf

Comments

Categories