Create custom field in Mechanical DPF

Options

Hello,

I have a script which I would like to use to plot certain parameter I calculate. To do that, I figured I could create a field and populate it in the loop. The script seems to work for simple geometry, please see picture below:

However, when I run it on actual geometry I get this

The script uses Property Provided to provide geometry selected by the user, in this case a single body. For some reason, when script is run on the more complicated geometry I get more bodies than 1 (originally the bottom plate is selected) and no results are displayed (they will not be 0).

Please see the script below:

def post_started(sender, analysis):# Do not edit this line
    define_dpf_workflow(analysis)

# Uncomment this function to enable retrieving results from the table/chart
# def table_retrieve_result(value):# Do not edit this line
    # import mech_dpf
    # import Ans.DataProcessing as dpf
    # wf = dpf.Workflow(this.WorkflowId)
    # wf.Connect('contour_selector', value)
    # this.Evaluate()

def define_dpf_workflow(analysis):
    import mech_dpf
    import Ans.DataProcessing as dpf
    mech_dpf.setExtAPI(ExtAPI)
    dataSource = dpf.DataSources(analysis.ResultFileName)
    meshData = DataModel.MeshDataByName("Global")
    timeScop = dpf.Scoping()
    timeScop.Ids = [3]
    timeScop.Location='Elemental'

    scoping_refs = this.GetCustomPropertyByPath("Geometry/Scoping Property/Geometry Selection").Value
    nodes_loc = mech_dpf.GetNodeScopingByRefId(scoping_refs)



    #s1= dpf.operators.result.stress_principal_1() 
    #s2= dpf.operators.result.stress_principal_3() 
    #s3= dpf.operators.result.stress_principal_3()

    s1= dpf.operators.result.stress_X() 
    s2= dpf.operators.result.stress_Y() 
    s3= dpf.operators.result.stress_Z() 

    seqv = dpf.operators.result.stress_von_mises()

    Add_1_2 = dpf.operators.math.add()
    Add_1_2_3 = dpf.operators.math.add()
    op_1_3 = dpf.operators.math.scale()
    op_div = dpf.operators.math.component_wise_divide()
    op_mult_scale = dpf.operators.math.scale()
    op_div_2 = dpf.operators.math.component_wise_divide()

    s1.inputs.time_scoping.Connect(timeScop)
    s2.inputs.time_scoping.Connect(timeScop)
    s3.inputs.time_scoping.Connect(timeScop)
    seqv.inputs.time_scoping.Connect(timeScop)

    s1.inputs.data_sources.Connect(dataSource)
    s2.inputs.data_sources.Connect(dataSource)
    s3.inputs.data_sources.Connect(dataSource)
    seqv.inputs.data_sources.Connect(dataSource)

    s1.inputs.mesh_scoping.Connect(nodes_loc)
    s2.inputs.mesh_scoping.Connect(nodes_loc)
    s3.inputs.mesh_scoping.Connect(nodes_loc)
    seqv.inputs.mesh_scoping.Connect(nodes_loc)

    Add_1_2.inputs.fieldA.Connect(s1)
    Add_1_2.inputs.fieldB.Connect(s2)
    s_1_2 = Add_1_2.outputs.field.GetData()

    Add_1_2_3.inputs.fieldA.Connect(s_1_2)
    Add_1_2_3.inputs.fieldB.Connect(s3)
    s_1_2_3 = Add_1_2_3.outputs.field.GetData()

    op_1_3.inputs.field.Connect(s_1_2_3)
    op_1_3.inputs.ponderation.Connect(0.3333333333)
    #op_1_3.inputs.boolean.Connect(my_boolean)# optional
    sigmaH = op_1_3.outputs.field.GetData()

    op_div.inputs.fieldA.Connect(sigmaH)
    op_div.inputs.fieldB.Connect(seqv)
    H = op_div.outputs.field.GetData()

    Min_tol_strain = 0.05
    Rm=470 #Mpa
    Rp = 345 #MPa
    hard_factor = 0.5*(1+Rm/Rp)
    A = 0.22
    E = 210000 #MPa
    ref_strain = A+Rp/E
    tol_strain = 0

    npl_field = dpf.FieldsFactory.CreateScalarField(H.Data.Count, "Nodal")

    input_support = dpf.FieldsFactory.Create3DVectorField(H.Data.Count, "Nodal")
    input_field = dpf.FieldsFactory.CreateScalarField(H.Data.Count, "Nodal")
    lineIndex=1



    for node_Id in H.Scoping.Ids:
        if H.GetEntityDataById(node_Id) <=1/3:
            tol_strain = ref_strain
        else: tol_strain = Min_tol_strain + ((ref_strain-Min_tol_strain)/3)^(3*H.GetEntityDataById(node_Id))
        node = meshData.NodeById(node_Id) 
        input_support.Add(lineIndex, [node.X, node.Y, node.Z])
        # npl_field.UpdateEntityDataByEntityIndex(i,MIN(sqrt(E*tol_strain/Rp),hard_factor))
        input_field.Add(lineIndex,[min(sqrt(E*tol_strain/Rp),hard_factor)])
        lineIndex += 1
        #npl_field.UpdateEntityDataByEntityIndex(i,10.00)

    #DPF Mapper Operator
    rbf_prep = dpf.Operator("make_rbf_mapper")
    rbf_prep.Connect(0, input_support)
    rbf_apply = dpf.Operator("apply_mapping")
    rbf_apply.Connect(0, rbf_prep)
    rbf_apply.Connect(1, input_field)
    rbf_apply.Connect(2, input_support)
    #Radius to read results for mapping
    r = 0.01
    rbf_prep.Connect(1, r)
    npl_field = rbf_apply.GetOutputAsField(0)


    op_mult_scale.inputs.field.Connect(npl_field)
    op_mult_scale.inputs.ponderation.Connect(Rp)
    sigmaSK = op_mult_scale.outputs.field.GetData()

    op_div_2.inputs.fieldA.Connect(seqv)
    op_div_2.inputs.fieldB.Connect(sigmaSK)
    UC = op_div_2.outputs.field.GetData()

    # timeScop = dpf.Scoping()
    # timeScop.Ids = [1]
    # u.inputs.time_scoping.Connect(timeScop)
    dpf_workflow = dpf.Workflow()
    dpf_workflow.Add(s1)
    dpf_workflow.Add(s2)
    dpf_workflow.Add(s3)
    dpf_workflow.Add(seqv)
    dpf_workflow.Add(Add_1_2)
    dpf_workflow.Add(Add_1_2_3)
    dpf_workflow.Add(op_1_3)
    dpf_workflow.Add(op_div)
    dpf_workflow.Add(op_mult_scale)
    dpf_workflow.Add(op_div_2)

    #dpf_workflow.Add(nrm)
    # dpf_workflow.SetInputName(u, 0, 'time')
    # dpf_workflow.Connect('time', timeScop)
    dpf_workflow.SetOutputContour(op_div_2)
    dpf_workflow.Record('wf_id', False)
    this.WorkflowId = dpf_workflow.GetRecordedId()

Operator op_div seems to display correctly so there is something going on with npl_field. To be honest, I am not sure if rbf mapper is the way to go but it makes sense in my head:
1. Create a field of nodes positions input_support (will be the same nodes as other H field's scoping
2. Loop over all the nodes and populate input_field
3. Outside of the loop, use rbf mapper to create npl_field. Output support is equal to input support (I am mapping over identical locations).
4. Feed npl_field into some operators to finish of the calculations.
5. Plot op_div_2 operator results.

Can you please help me figure out if this is the correct approach and indicate why this script doesn't work for larger geometry? Larger geometry also has contacts and preloaded bolts so I wonder if this has some effect? Unfortunately I cannot share geometry model.

Answers