.. DO NOT EDIT. .. THIS FILE WAS AUTOMATICALLY GENERATED BY SPHINX-GALLERY. .. TO MAKE CHANGES, EDIT THE SOURCE PYTHON FILE: .. "examples\09-averaging\01-average_across_bodies.py" .. LINE NUMBERS ARE GIVEN BELOW. .. only:: html .. note:: :class: sphx-glr-download-link-note :ref:`Go to the end ` to download the full example code. .. rst-class:: sphx-glr-example-title .. _sphx_glr_examples_09-averaging_01-average_across_bodies.py: .. _ref_average_across_bodies: Average across bodies ~~~~~~~~~~~~~~~~~~~~~ In multibody simulations, some nodes may be shared by the bodies at their interfaces, but the values of the results (for example stresses or strains) calculated at these nodes may differ between the bodies. This can cause discontinuous plots, given that a single node will have multiple values for a variable. To avoid this, you can average these results across the bodies of the model. This example demonstrates how to average across bodies in DPF when dealing with ``Nodal`` variables. It also illustrates how the end results of a postprocessing workflow can be different when averaging and when not. .. note:: This example requires DPF 6.1 or above. For more information, see :ref:`ref_compatibility`. .. GENERATED FROM PYTHON SOURCE LINES 46-47 Import the necessary modules .. GENERATED FROM PYTHON SOURCE LINES 47-53 .. code-block:: Python from ansys.dpf import core as dpf from ansys.dpf.core import operators as ops from ansys.dpf.core import examples .. GENERATED FROM PYTHON SOURCE LINES 54-55 Load the simulation results from an RST file and create a model of it. .. GENERATED FROM PYTHON SOURCE LINES 55-60 .. code-block:: Python analysis = examples.download_piston_rod() model = dpf.Model(analysis) print(model) .. rst-class:: sphx-glr-script-out .. code-block:: none DPF Model ------------------------------ Static analysis Unit system: MKS: m, kg, N, s, V, A, degC Physics Type: Mechanical Available results: - displacement: Nodal Displacement - reaction_force: Nodal Force - stress: ElementalNodal Stress - elemental_volume: Elemental Volume - stiffness_matrix_energy: Elemental Energy-stiffness matrix - artificial_hourglass_energy: Elemental Hourglass Energy - thermal_dissipation_energy: Elemental thermal dissipation energy - kinetic_energy: Elemental Kinetic Energy - co_energy: Elemental co-energy - incremental_energy: Elemental incremental energy - elastic_strain: ElementalNodal Strain - element_euler_angles: ElementalNodal Element Euler Angles - structural_temperature: ElementalNodal Structural temperature ------------------------------ DPF Meshed Region: 33337 nodes 18235 elements Unit: m With solid (3D) elements ------------------------------ DPF Time/Freq Support: Number of sets: 3 Cumulative Time (s) LoadStep Substep 1 1.000000 1 1 2 2.000000 2 1 3 3.000000 3 1 .. GENERATED FROM PYTHON SOURCE LINES 61-64 To visualize the model and see how the bodies are connected, extract their individual meshes using the ``split_mesh`` operator with the ``mat`` (or "material") property. .. GENERATED FROM PYTHON SOURCE LINES 64-71 .. code-block:: Python mesh = model.metadata.meshed_region split_mesh_op = ops.mesh.split_mesh(mesh=mesh, property="mat") meshes = split_mesh_op.outputs.meshes() meshes.plot(text="Body meshes") .. image-sg:: /examples/09-averaging/images/sphx_glr_01-average_across_bodies_001.png :alt: 01 average across bodies :srcset: /examples/09-averaging/images/sphx_glr_01-average_across_bodies_001.png :class: sphx-glr-single-img .. GENERATED FROM PYTHON SOURCE LINES 72-74 As can be seen in the preceding image, even though the piston rod is one single part, it is composed of two different bodies. Additionally, their interface shares common nodes. .. GENERATED FROM PYTHON SOURCE LINES 76-83 Averaging across bodies with DPF ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ # To compare the results of averaging across bodies and not averaging, define two workflows. The variable of interest is the Von Mises stress field, which is calculated by applying the ``eqv_fc`` operator on the stresses extracted from the model. .. GENERATED FROM PYTHON SOURCE LINES 85-119 .. graphviz:: digraph foo { graph [pad="0", nodesep="0.3", ranksep="0.3"] node [shape=box, style=filled, fillcolor="#ffcc0", margin="0"]; rankdir=LR; splines=line; node [fixedsize=true,width=2.5] ds [label="data_src", shape=box, style=filled, fillcolor=cadetblue2]; stress [label="stress"]; scp [label="split_on_property_type"]; eln_to_n ["elemental_nodal_to_nodal_fc"]; vm [label="eqv_fc"]; avg [label="weighted_merge_fields_by_label"]; subgraph cluster_1 { ds -> scp [style=dashed]; scp -> stress; stress -> eln_to_n; eln_to_n -> vm; label="Without averaging across bodies"; style=filled; fillcolor=lightgrey; } subgraph cluster_2 { ds -> scp [style=dashed]; scp -> stress; stress -> eln_to_n; eln_to_n -> vm; vm -> avg; label="With averaging across bodies"; style=filled; fillcolor=lightgrey; } } .. GENERATED FROM PYTHON SOURCE LINES 122-128 Workflow for not averaging across bodies ---------------------------------------- Computing Von Mises stresses without averaging across the bodies of the model requires the stresses to be extracted separately for each body. To do this in DPF, pass a scopings container the stress operator that contains the elements of each body in scopings, separated by the ``mat`` label .. GENERATED FROM PYTHON SOURCE LINES 128-136 .. code-block:: Python split_scop_op = ops.scoping.split_on_property_type() split_scop_op.inputs.mesh.connect(mesh) split_scop_op.inputs.requested_location.connect(dpf.locations.elemental) split_scop_op.inputs.label1.connect("mat") print(split_scop_op.outputs.mesh_scoping()) .. rst-class:: sphx-glr-script-out .. code-block:: none DPF Scopings Container with 2 scoping(s) defined on labels: mat with: - scoping 0 {mat: 1, } located on Elemental 9026 entities. - scoping 1 {mat: 2, } located on Elemental 9209 entities. .. GENERATED FROM PYTHON SOURCE LINES 137-138 Set the time set of interest to the last time set: .. GENERATED FROM PYTHON SOURCE LINES 138-140 .. code-block:: Python time_set = 3 .. GENERATED FROM PYTHON SOURCE LINES 141-142 Extracting the stresses for each body of the simulation: .. GENERATED FROM PYTHON SOURCE LINES 142-149 .. code-block:: Python stress_op = ops.result.stress() stress_op.inputs.time_scoping.connect(time_set) stress_op.inputs.data_sources.connect(model) stress_op.inputs.mesh_scoping.connect(split_scop_op) stress_op.inputs.requested_location.connect(dpf.locations.elemental_nodal) .. GENERATED FROM PYTHON SOURCE LINES 150-151 Proceeding with the workflow to obtain ``Nodal`` Von Mises stresses: .. GENERATED FROM PYTHON SOURCE LINES 151-159 .. code-block:: Python eln_to_n_op = ops.averaging.elemental_nodal_to_nodal_fc() eln_to_n_op.inputs.fields_container.connect(stress_op) von_mises_op = ops.invariant.von_mises_eqv_fc() von_mises_op.inputs.fields_container.connect(eln_to_n_op) print(von_mises_op.outputs.fields_container()) .. rst-class:: sphx-glr-script-out .. code-block:: none DPF stress(s)Fields Container with 2 field(s) defined on labels: mat time with: - field 0 {mat: 1, time: 3} with Nodal location, 1 components and 17281 entities. - field 1 {mat: 2, time: 3} with Nodal location, 1 components and 17610 entities. .. GENERATED FROM PYTHON SOURCE LINES 160-164 As you can see, the final Von Mises stresses fields container has the ``mat`` label with two different entries, meaning that it holds data for two separate bodies. Finally, define this workflow as a function for better organization and ease of use: .. GENERATED FROM PYTHON SOURCE LINES 164-197 .. code-block:: Python def not_average_across_bodies(analysis): # This function extracts the ElementalNodal stress tensors of the simulation # for each body involved, averages them to the nodes and computes Von Mises model = dpf.Model(analysis) mesh = model.metadata.meshed_region time_set = 3 split_scop_op = ops.scoping.split_on_property_type() split_scop_op.inputs.mesh.connect(mesh) split_scop_op.inputs.requested_location.connect(dpf.locations.elemental) split_scop_op.inputs.label1.connect("mat") stress_op = ops.result.stress() stress_op.inputs.time_scoping.connect(time_set) stress_op.inputs.data_sources.connect(model) stress_op.inputs.mesh_scoping.connect(split_scop_op) stress_op.inputs.requested_location.connect(dpf.locations.elemental_nodal) eln_to_n_op = ops.averaging.elemental_nodal_to_nodal_fc() eln_to_n_op.inputs.fields_container.connect(stress_op) von_mises_op = ops.invariant.von_mises_eqv_fc() von_mises_op.inputs.fields_container.connect(eln_to_n_op) vm_stresses = von_mises_op.outputs.fields_container() return vm_stresses .. GENERATED FROM PYTHON SOURCE LINES 198-205 Workflow for averaging across bodies ------------------------------------ The workflow for performing averaging across bodies in DPF is similar to to the one shown above, with the extraction of stresses per body. The difference comes in the end, where a weighted merge is done between the fields that contain different values for the ``mat`` label to actually average the results across the bodies. Define a function like the one above: .. GENERATED FROM PYTHON SOURCE LINES 205-247 .. code-block:: Python def average_across_bodies(analysis): # This function extracts the ElementalNodal stress tensors of the simulation # for each body involved, averages them to the nodes and computes Von Mises model = dpf.Model(analysis) mesh = model.metadata.meshed_region time_set = 3 split_scop_op = ops.scoping.split_on_property_type() split_scop_op.inputs.mesh.connect(mesh) split_scop_op.inputs.requested_location.connect(dpf.locations.elemental) split_scop_op.inputs.label1.connect("mat") stress_op = ops.result.stress() stress_op.inputs.time_scoping.connect(time_set) stress_op.inputs.data_sources.connect(model) stress_op.inputs.mesh_scoping.connect(split_scop_op) stress_op.inputs.requested_location.connect(dpf.locations.elemental_nodal) eln_to_n_op = ops.averaging.elemental_nodal_to_nodal_fc() eln_to_n_op.inputs.fields_container.connect(stress_op) # Mid node weights needed for averaging across bodies eln_to_n_op.inputs.extend_weights_to_mid_nodes.connect(True) von_mises_op = ops.invariant.von_mises_eqv_fc() von_mises_op.inputs.fields_container.connect(eln_to_n_op) # Merging fields that represent different bodies merge_op = ops.utility.weighted_merge_fields_by_label() merge_op.inputs.fields_container.connect(von_mises_op) merge_op.inputs.label.connect("mat") # Connecting weights needed to perform the weighted average merge_op.connect(1000, eln_to_n_op, 1) vm_stresses = merge_op.outputs.fields_container() return vm_stresses .. GENERATED FROM PYTHON SOURCE LINES 248-250 In this case, we can see that the output fields container only has one field, indicating that the results of the two different bodies were averaged successfully. .. GENERATED FROM PYTHON SOURCE LINES 250-253 .. code-block:: Python print(average_across_bodies(analysis)) .. rst-class:: sphx-glr-script-out .. code-block:: none DPF Fields Container with 1 field(s) defined on labels: time with: - field 0 {time: 3} with Nodal location, 1 components and 33337 entities. .. GENERATED FROM PYTHON SOURCE LINES 254-259 Plot and compare the results ~~~~~~~~~~~~~~~~~~~~~~~~~~~~ The two different approaches can be compared. The first plot shows the results when averaging across bodies is not performed, while the second illustrates when it is. .. GENERATED FROM PYTHON SOURCE LINES 259-265 .. code-block:: Python non_avg_stresses = not_average_across_bodies(analysis) avg_stresses = average_across_bodies(analysis) meshes.plot(non_avg_stresses) mesh.plot(avg_stresses) .. rst-class:: sphx-glr-horizontal * .. image-sg:: /examples/09-averaging/images/sphx_glr_01-average_across_bodies_002.png :alt: 01 average across bodies :srcset: /examples/09-averaging/images/sphx_glr_01-average_across_bodies_002.png :class: sphx-glr-multi-img * .. image-sg:: /examples/09-averaging/images/sphx_glr_01-average_across_bodies_003.png :alt: 01 average across bodies :srcset: /examples/09-averaging/images/sphx_glr_01-average_across_bodies_003.png :class: sphx-glr-multi-img .. GENERATED FROM PYTHON SOURCE LINES 266-267 Finally, the maximum stresses for both cases can be compared: .. GENERATED FROM PYTHON SOURCE LINES 267-287 .. code-block:: Python min_max = dpf.operators.min_max.min_max_fc() # Non averaged across bodies min_max.inputs.fields_container.connect(non_avg_stresses) max_non_avg = max(min_max.outputs.field_max().data) # Averaged across bodies min_max.inputs.fields_container.connect(avg_stresses) max_avg = max(min_max.outputs.field_max().data) diff = abs(max_avg - max_non_avg) / max_non_avg * 100 print("Max stress when averaging across bodies is activated: {:.2f} Pa".format(max_avg)) print("Max stress when averaging across bodies is deactivated: {:.2f} Pa".format(max_non_avg)) print( "The maximum stress value when averaging across bodies is PERFORMED \ is {:.2f}% LOWER than when it is NOT PERFORMED".format(diff) ) .. rst-class:: sphx-glr-script-out .. code-block:: none Max stress when averaging across bodies is activated: 10786714384.57 Pa Max stress when averaging across bodies is deactivated: 10796069601.08 Pa The maximum stress value when averaging across bodies is PERFORMED is 0.09% LOWER than when it is NOT PERFORMED .. GENERATED FROM PYTHON SOURCE LINES 288-298 Dedicated Operator ~~~~~~~~~~~~~~~~~~ .. note:: The operator detailed below is available in Ansys 23R2 and later versions. Alternatively, those workflows can be automatically instantiated by calling the ``stress_eqv_as_mechanical`` operator, which does exactly the same thing as described in the functions above, depending on what is passed to the "average_across_bodies" input pin: .. GENERATED FROM PYTHON SOURCE LINES 298-306 .. code-block:: Python stress_op = ops.result.stress_eqv_as_mechanical() stress_op.inputs.time_scoping.connect([time_set]) stress_op.inputs.data_sources.connect(model) stress_op.inputs.requested_location.connect(dpf.locations.nodal) stress_op.inputs.average_across_bodies.connect(False) print(stress_op.outputs.fields_container()) .. rst-class:: sphx-glr-script-out .. code-block:: none DPF stress(s)Fields Container with 2 field(s) defined on labels: mat time with: - field 0 {mat: 1, time: 3} with Nodal location, 1 components and 17281 entities. - field 1 {mat: 2, time: 3} with Nodal location, 1 components and 17610 entities. .. rst-class:: sphx-glr-timing **Total running time of the script:** (0 minutes 15.569 seconds) .. _sphx_glr_download_examples_09-averaging_01-average_across_bodies.py: .. only:: html .. container:: sphx-glr-footer sphx-glr-footer-example .. container:: sphx-glr-download sphx-glr-download-jupyter :download:`Download Jupyter notebook: 01-average_across_bodies.ipynb <01-average_across_bodies.ipynb>` .. container:: sphx-glr-download sphx-glr-download-python :download:`Download Python source code: 01-average_across_bodies.py <01-average_across_bodies.py>` .. container:: sphx-glr-download sphx-glr-download-zip :download:`Download zipped: 01-average_across_bodies.zip <01-average_across_bodies.zip>` .. only:: html .. rst-class:: sphx-glr-signature `Gallery generated by Sphinx-Gallery `_