Python Result fatigue calculation

Pernelle Marone-Hitz

Is there an example of fatigue calculation with DPF?

Pernelle Marone-Hitz

The following example has no physical meaning and is just an illustration of DPF capabilies.
It uses the dummy custom stress computed in this example, and reads the S-N curve of the material of the first body in the tree (adapt to your use case!).

Note - a value of -1 or -2 in the result plot will indicate that the interpolation of the number of cycles to failure could not be performed (usually due to lack of material data to interpolate).

Complete code:

def LinearInterpolation(x,x1,y1,x2,y2):
    # Simple linear interpolation
    # Returns -1 if x is outside of [x1,x2]
    if x<x1 or x>x2:
        return -1 # out of range
    if y1==y2:
        return y1
    else:
        a=(y2-y1)/(x2-x1)
        b=y1-a*x1
        return a*x+b

def GetLocalCyclesToFailure(sn_curve, local_stress_si):
    cycles_list = sn_curve["Cycles"][1:]
    alternating_stress_list = sn_curve["Alternating Stress"][1:]
    index = 0
    # Find closest cycle number in stress_data
    for i in range(len(cycles_list)):
        alt_stress = alternating_stress_list[i]
        if alt_stress < local_stress_si:
            index=i
            break
    # Get (cycles, alternating stress) points between which to interpolate
    if index!= len(cycles_list):
        y1=cycles_list[index]
        y2=cycles_list[index-1]
        x1=alternating_stress_list[index]
        x2=alternating_stress_list[index-1]
    else:
        return -2 # value is out of bounds
    return LinearInterpolation(local_stress_si,x1,y1,x2,y2) # return number of cycles

def GetFieldCyclesToFailure(stress_op, sn_curve, analysis):
    import mech_dpf
    import Ans.DataProcessing as dpf
    node_ids = stress_op.outputs.getfields_container()[0].ScopingIds
    limit_field = dpf.FieldsFactory.CreateScalarField(1,'Nodal')
    stress_data = stress_op.outputs.getfields_container()[0]
    for nid in node_ids:
        local_stress =stress_data.GetEntityDataById(nid)[0]
        # Convert local stress to SI units
        local_stress_quantity = Quantity(str(local_stress) + '[' + analysis.CurrentConsistentUnitFromQuantityName("Stress") + ']')
        local_stress_si = local_stress_quantity.ConvertToUnitSystem("SI", "Stress").Value
        limit_field.Add(nid,[float(GetLocalCyclesToFailure(sn_curve, local_stress_si))])
    return limit_field


def define_dpf_workflow(analysis):
    import mech_dpf
    import Ans.DataProcessing as dpf
    mech_dpf.setExtAPI(ExtAPI)
    data_source = dpf.DataSources(analysis.ResultFileName)
    model = dpf.Model(data_source)
    mesh =  model.Mesh

    # Get SN curves for material defined on first body in the Mechanical tree
    import materials
    first_body = Model.GetChildren(DataModelObjectCategory.Body, True)[0]
    material = ExtAPI.DataModel.GetObjectsByName(first_body.Material)[0]
    mat_enginerring_data = material.GetEngineeringDataMaterial()
    sn_curve = materials.GetMaterialPropertyByName(mat_enginerring_data,"S-N Curve")

    # Time sets
    timePointsOp =dpf.operators.metadata.time_freq_provider()
    timePointsOp.inputs.data_sources.Connect(data_source)
    timepoints = dpf.operators.metadata.time_freq_support_get_attribute(time_freq_support=timePointsOp,property_name="time_freqs",)

    # S1
    prin1 = dpf.operators.result.stress_principal_1()
    prin1.inputs.data_sources.Connect(data_source)
    prin1.inputs.time_scoping.Connect(timepoints)
    # S3
    prin3 = dpf.operators.result.stress_principal_3()
    prin3.inputs.data_sources.Connect(data_source)
    prin3.inputs.time_scoping.Connect(timepoints)
    # S1-S3
    minus_fc = dpf.operators.math.minus_fc()
    minus_fc.inputs.field_or_fields_container_A.Connect(prin1.outputs.fields_container)
    minus_fc.inputs.field_or_fields_container_B.Connect(prin3.outputs.fields_container)
    # Scale factor
    scale = dpf.operators.math.scale_fc()
    scale.inputs.fields_container.Connect(minus_fc)
    scale.inputs.ponderation.Connect(sqrt(2)/2)

    # Get alternate stress
    min_max = dpf.operators.min_max.min_max_over_time_by_entity()
    min_max.inputs.fields_container.Connect(scale)
    min_field_by_time =  min_max.outputs.getmin()
    max_field_by_time =  min_max.outputs.getmax()
    alternate_stress = dpf.operators.math.minus_fc()
    alternate_stress.inputs.field_or_fields_container_A.Connect(min_max.outputs.getmax())
    alternate_stress.inputs.field_or_fields_container_B.Connect(min_max.outputs.getmin())

    # Get cycles to failure
    new_field = GetFieldCyclesToFailure(alternate_stress, sn_curve, analysis)
    output = dpf.operators.utility.forward()
    output.inputs.any.Connect(new_field)

    # Plot
    dpf_workflow = dpf.Workflow()
    dpf_workflow.Add(output)
    dpf_workflow.SetOutputContour(output)
    dpf_workflow.Record('wf_id', False)
    this.WorkflowId = dpf_workflow.GetRecordedId()