How do I calculate linearized stress (Membrane, Bending and Peak) using DPF?

Options
Ayush Kumar
Ayush Kumar Member, Moderator, Employee Posts: 414
First Anniversary Ansys Employee Solution Developer Community of Practice Member First Answer
edited June 2023 in Structures

How do I calculate linearized stress (Membrane, Bending and Peak) using DPF?

Tagged:

Comments

  • Ayush Kumar
    Ayush Kumar Member, Moderator, Employee Posts: 414
    First Anniversary Ansys Employee Solution Developer Community of Practice Member First Answer
    edited May 2023
    Options
    import mech_dpf
    import Ans.DataProcessing as dpf
    
    
    path = ExtAPI.DataModel.GetObjectsByName("Path")[0]
    x_0 = path.StartXCoordinate.Value
    y_0 = path.StartYCoordinate.Value
    z_0 = path.StartZCoordinate.Value
    
    
    x_1 = path.EndXCoordinate.Value
    y_1 = path.EndYCoordinate.Value
    z_1 = path.EndZCoordinate.Value
    
    
    path_length = ((x_1 - x_0) ** 2 + (y_1 - y_0) ** 2 + (z_1 - z_0) ** 2) ** 0.5
    
    
    n_points = path.DiscretizationPoints + 2
    delta = path_length / (n_points - 1)
    
    
    line_unit_vector = [(x_1 - x_0) / path_length, (y_1 - y_0) / path_length, (z_1 - z_0) / path_length]
    
    
    # Line equation
    fx = lambda t: x_0 + line_unit_vector[0] * t
    fy = lambda t: y_0 + line_unit_vector[1] * t
    fz = lambda t: z_0 + line_unit_vector[2] * t
    
    
    coordinates = [[fx(t * delta), fy(t * delta), fz(t * delta)] for t in range(n_points)]
    flat_coordinates = [entry for data in coordinates for entry in data]
    
    
    field_coord = dpf.FieldsFactory.Create3DVectorField(n_points)
    field_coord.Data = flat_coordinates
    field_coord.Scoping.Ids = list(range(1, n_points + 1))
    
    
    model = dpf.Model(Model.Analyses[0].ResultFileName)
    
    
    ts = dpf.Scoping()
    ts.Ids = [2]
    
    
    # SX
    s = model.CreateOperator("SX")
    s.inputs.requested_location.Connect("Nodal")
    s.inputs.time_scoping.Connect(ts)
    s2_f = s.outputs.fields_container.GetData()
    
    
    mapping_operator = dpf.operators.mapping.on_coordinates(
        fields_container=s2_f,
        coordinates=field_coord,
        create_support=True,
        mesh=model.Mesh
    )
    fields_mapped = mapping_operator.outputs.fields_container.GetData()
    
    
    ls = list(fields_mapped[0].Data)
    # Membrane Stress
    membrane_stress = (ls[0] / 2 + ls[-1] / 2 + sum(ls[1:-1])) / (n_points - 1)
    
    
    # Bending stress
    path_1 = -1 * path_length / 2
    path_n = path_length / 2
    path_range = [path_1 + delta * i for i in range(n_points)]
    path_range_field = dpf.FieldsFactory.CreateScalarField(numEntities=n_points, location=dpf.locations.nodal)
    path_range_field.Data = path_range
    path_range_field.ScopingIds = range(1, 50)
    
    
    # Function to be integrated 
    stress_scaled = dpf.operators.math.scale_by_field(fieldA=fields_mapped[0], fieldB=path_range_field)
    stress_scaled_data = list(stress_scaled.outputs.field.GetData().Data)
    
    
    # Use extended Simpson's rule for Numerical Integration of `stress_scaled_data`
    stress_scaled_integral = (17 * stress_scaled_data[0] + 59 * stress_scaled_data[1] + 43 * stress_scaled_data[2] + 49 *
                              stress_scaled_data[3] + 48 * sum(stress_scaled_data[4:-4]) + 49 * stress_scaled_data[
                                  n_points - 4] + 43 * stress_scaled_data[n_points - 3] + 59 * stress_scaled_data[
                                  n_points - 2] + 17 * stress_scaled_data[n_points - 1]) / 48.0
    
    
    b1 = stress_scaled_integral * (-6.0 / path_length) / 48.0
    b2 = (-1.0) * b1
    <?WORD FLAGGED?>_bending = (b2 - b1) / (path_length)
    fb = lambda t: b1 + <?WORD FLAGGED?>_bending * t
    bending_stress = [fb(t * delta) for t in range(n_points)]
    
    
    bend_f = dpf.FieldsFactory.CreateScalarField(numEntities=49, location=dpf.locations.nodal)
    bend_f.Data = bending_stress
    bend_f.Scoping.Ids = list(range(1, n_points + 1))
    
    
    # Add Membrane and Bending Stress
    bend_mem_f = dpf.operators.math.add_constant(field=bend_f, ponderation=membrane_stress)
    bend_mem_f_neg = dpf.operators.math.scale(field=bend_mem_f, ponderation=-1.0)
    
    
    # Peak stress
    ps = dpf.operators.math.add(fieldA=fields_mapped[0], fieldB=bend_mem_f_neg)
    
    
    
    
    
  • dafedin
    dafedin Member Posts: 21
    First Anniversary First Comment Name Dropper
    Options

    Hello Ayush,

    Thanks for this comprehensive example. I wonder if there is a way to solve the reversed task: find elements of the model whose results were used for mapping to path points. Maybe there is another option rather than screening over the whole mesh and finding elements whose centroids are situated in the path direction.

    Regards

  • Ayush Kumar
    Ayush Kumar Member, Moderator, Employee Posts: 414
    First Anniversary Ansys Employee Solution Developer Community of Practice Member First Answer
    Options

    Hi @dafedin,

    That information is not provided as an output of the operator at the moment.