Create custom field in Mechanical DPF

ArnoldH · June 19

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.

Rohith Patchigolla · July 15

Hello @ArnoldH

Please find a modified script here for your reference.
Some comments on changes:
• When defining integer numbers which may be involved in a division later on, always add .0 after the number (Eg: Rm = 470.0), else the division will lead to rounded integer
• Time scoping corresponds to result set number and not the time.
• rbf map operator is not needed
• To plot a field created in a python result, please use forward_field operator
• For scoping of the new field, I used the scoping IDs coming in directly from H to maintain the correct order of nodes.

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
    import math

    mech_dpf.setExtAPI(ExtAPI)
    dataSource = dpf.DataSources(analysis.ResultFileName)

    timeScop = dpf.Scoping()
    timeScop.Ids = [3]

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

    #ns_op = dpf.operators.scoping.on_named_selection(data_sources = dataSource, requested_location = "Nodal", named_selection_name = "SELECTION")

    s1= dpf.operators.result.stress_X(data_sources = dataSource, time_scoping = timeScop, mesh_scoping = nodes_loc) 
    s2= dpf.operators.result.stress_Y(data_sources = dataSource, time_scoping = timeScop, mesh_scoping = nodes_loc) 
    s3= dpf.operators.result.stress_Z(data_sources = dataSource, time_scoping = timeScop, mesh_scoping = nodes_loc) 
    seqv = dpf.operators.result.stress_von_mises(data_sources = dataSource, time_scoping = timeScop, mesh_scoping = nodes_loc)

    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()
    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)
    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.0 #Mpa
    Rp = 345.0 #MPa
    hard_factor = 0.5*(1+Rm/Rp)
    A = 0.22
    E = 210000 #MPa
    ref_strain = A+Rp/E
    tol_strain = 1

    npl_field = dpf.FieldsFactory.CreateScalarField(H.Data.Count, "Nodal")
    npl_field_Data = []
    for node_Id in H.Scoping.Ids:
        if H.GetEntityDataById(node_Id)[0] <=1/3:
            tol_strain = ref_strain
        else:
            tol_strain = Min_tol_strain + math.pow((ref_strain-Min_tol_strain)/3,(3*H.GetEntityDataById(node_Id)[0]))
        npl_field_Data.append(min(sqrt(E*tol_strain/Rp),hard_factor))

    npl_field.Data = npl_field_Data
    #npl_field.Data = H.Data

    my_scoping = dpf.Scoping()
    my_scoping.Location = "Nodal"
    my_scoping.Ids = H.Scoping.Ids
    npl_field.Scoping = my_scoping

    # forward field to an operator to plot it
    forward = dpf.operators.utility.forward_field()
    forward.inputs.field.Connect(npl_field)

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

Chris Harrold · June 24

Hello, @ArnoldH - this post looks to be (mostly) a duplicate of this post:
https://discuss.ansys.com/discussion/3410

If that is right, can I remove this one as a duplicate and stick with the other one, or vice/versa?

Thanks!
Chris

ArnoldH · June 25

Hi @Chris Harrold , sorry for the duplicate. Can you remove the other post, I edited this one a bit more. Thanks

Chris Harrold · June 26

Will do, @ArnoldH and @AKD-Scripting-Team might be able to give you some guidance here.

Create custom field in Mechanical DPF

Best Answer

Answers

Categories