.. DO NOT EDIT. .. THIS FILE WAS AUTOMATICALLY GENERATED BY SPHINX-GALLERY. .. TO MAKE CHANGES, EDIT THE SOURCE PYTHON FILE: .. "auto_examples/plates/plot_plate08_triangle_patch_test.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_auto_examples_plates_plot_plate08_triangle_patch_test.py: ========================================================== Patch test for triangular plate under in-plane loading ========================================================== The patch test is an empirical minimum test which every finite element has to pass to ensure convergence with mesh refinement. It consists of a problem for which a known homogeneous solution exists. For plates, we commonly use a rectangular plate subject to homogeneous edge loading, e.g., constant tension in the x-direction, or constant shear, etc. The mesh must contain distorted elements and at least one element not attached to any node on the boundary. Author: Peter Mackenzie-Helnwein .. GENERATED FROM PYTHON SOURCE LINES 16-181 .. code-block:: Python import numpy as np from femedu.examples import Example from femedu.domain import System, Node from femedu.solver import NewtonRaphsonSolver from femedu.elements.linear import Triangle from femedu.materials import PlaneStress class ExamplePlate03(Example): def problem(self): # ========== setting mesh parameters ============== N = 8 # number of elements in the mesh Lx = 100.0 # length of plate in the x-direction Ly = 80.0 # length of plate in the y-direction # ========== setting material parameters ============== params = dict( E = 20000., # Young's modulus nu = 0.250, # Poisson's ratio t = 1.00 # thickness of the plate ) # ========== setting load parameters ============== px = 10.0 # uniform load normal to x=const py = 0.0 # uniform load normal to y=const pxy = 0.0 # uniform shear load on x=const and y=const # ========== setting analysis parameters ============== target_load_level = 1.00 # reference load max_steps = 2 # number of load steps: 2 -> [0.0, 1.0] # define a list of target load levels load_levels = np.linspace(0, target_load_level, max_steps) # # ==== Build the system model ==== # model = System() model.setSolver(NewtonRaphsonSolver()) # create nodes nodes = ( Node(0.0*Lx, 0.0*Ly), # nd 0 Node(0.2*Lx, 0.0*Ly), # nd 1 Node(0.5*Lx, 0.0*Ly), # nd 2 Node(0.7*Lx, 0.0*Ly), # nd 3 Node(1.0*Lx, 0.0*Ly), # nd 4 # Node(0.0*Lx, 0.2*Ly), # nd 5 Node(0.15*Lx,0.3*Ly), # nd 6 Node(0.5*Lx, 0.2*Ly), # nd 7 Node(0.8*Lx, 0.3*Ly), # nd 8 Node(1.0*Lx, 0.2*Ly), # nd 9 # Node(0.0*Lx, 0.6*Ly), # nd 10 Node(0.2*Lx, 0.5*Ly), # nd 11 Node(0.7*Lx, 0.7*Ly), # nd 12 Node(0.9*Lx, 0.6*Ly), # nd 13 Node(1.0*Lx, 0.8*Ly), # nd 14 # Node(0.0*Lx, 1.0*Ly), # nd 15 Node(0.3*Lx, 1.0*Ly), # nd 16 Node(0.55*Lx,1.0*Ly), # nd 17 Node(0.8*Lx, 1.0*Ly), # nd 18 Node(1.0*Lx, 1.0*Ly), # nd 19 ) elements = ( Triangle(nodes[0],nodes[1],nodes[5],PlaneStress(params)), # elem 0 Triangle(nodes[1],nodes[2],nodes[6],PlaneStress(params)), # elem 1 Triangle(nodes[2],nodes[3],nodes[7],PlaneStress(params)), # elem 2 Triangle(nodes[3],nodes[4],nodes[8],PlaneStress(params)), # elem 3 # Triangle(nodes[6],nodes[5],nodes[1],PlaneStress(params)), # elem 4 Triangle(nodes[7],nodes[6],nodes[2],PlaneStress(params)), # elem 5 Triangle(nodes[8],nodes[7],nodes[3],PlaneStress(params)), # elem 6 Triangle(nodes[9],nodes[8],nodes[4],PlaneStress(params)), # elem 7 # Triangle(nodes[5],nodes[6],nodes[10],PlaneStress(params)), # elem 8 Triangle(nodes[6],nodes[7],nodes[11],PlaneStress(params)), # elem 9 Triangle(nodes[7],nodes[8],nodes[12],PlaneStress(params)), # elem 10 Triangle(nodes[8],nodes[9],nodes[13],PlaneStress(params)), # elem 11 # Triangle(nodes[11],nodes[10],nodes[6],PlaneStress(params)), # elem 12 Triangle(nodes[12],nodes[11],nodes[7],PlaneStress(params)), # elem 13 Triangle(nodes[13],nodes[12],nodes[8],PlaneStress(params)), # elem 14 Triangle(nodes[14],nodes[13],nodes[9],PlaneStress(params)), # elem 15 # Triangle(nodes[10],nodes[11],nodes[15],PlaneStress(params)), # elem 16 Triangle(nodes[11],nodes[12],nodes[16],PlaneStress(params)), # elem 17 Triangle(nodes[12],nodes[13],nodes[17],PlaneStress(params)), # elem 18 Triangle(nodes[13],nodes[14],nodes[18],PlaneStress(params)), # elem 19 # Triangle(nodes[16],nodes[15],nodes[11],PlaneStress(params)), # elem 20 Triangle(nodes[17],nodes[16],nodes[12],PlaneStress(params)), # elem 21 Triangle(nodes[18],nodes[17],nodes[13],PlaneStress(params)), # elem 22 Triangle(nodes[19],nodes[18],nodes[14],PlaneStress(params)), # elem 23 ) model.addNode(*nodes) model.addElement(*elements) # define support(s) fix_x = (0,) fix_y = (0,4) for idx in fix_x: nodes[idx].fixDOF('ux') # horizontal support left end for idx in fix_y: nodes[idx].fixDOF('uy') # vertical support right end # ==== complete the reference load ==== # surface loads on the left side elements[ 0].setSurfaceLoad(2,px) elements[ 8].setSurfaceLoad(2,px) elements[16].setSurfaceLoad(2,px) # surface loads on the right side elements[ 7].setSurfaceLoad(2,px) elements[15].setSurfaceLoad(2,px) elements[23].setSurfaceLoad(2,px) # these are only nodal forces as part of the reference load # .. load only the upper node #print(model) model.report() model.plotSystem() model.setLoadFactor(10.0) model.solve() model.solver.showKt() model.report() model.plot(factor=25., filename="plate08_deformed.png") # tmpl = "{:12.4f} {:12.4f}" # for elem in model.elements: # print(tmpl.format(elem.area, elem.material.getStress()['xx'])) # requires femedu-1.0.25 or newer model.valuePlot('sxx', show_mesh=True) # tmpl = "{:>12} {:>5} {:12.4f} {:12.4f} {:12.4f}" # for node in model.nodes: # print(tmpl.format(node.getID(), node._mapped_variable, node._weighted_value, node._weight, node.getMappedValue('sxx'))) model.valuePlot('syy', show_mesh=True) model.valuePlot('sxy', show_mesh=True) .. GENERATED FROM PYTHON SOURCE LINES 199-201 Run the example by creating an instance of the problem and executing it by calling :py:meth:`Example.run()` .. GENERATED FROM PYTHON SOURCE LINES 201-207 .. code-block:: Python if __name__ == "__main__": ex = ExamplePlate03() ex.run() .. rst-class:: sphx-glr-horizontal * .. image-sg:: /auto_examples/plates/images/sphx_glr_plot_plate08_triangle_patch_test_001.png :alt: Undeformed system :srcset: /auto_examples/plates/images/sphx_glr_plot_plate08_triangle_patch_test_001.png :class: sphx-glr-multi-img * .. image-sg:: /auto_examples/plates/images/sphx_glr_plot_plate08_triangle_patch_test_002.png :alt: plot plate08 triangle patch test :srcset: /auto_examples/plates/images/sphx_glr_plot_plate08_triangle_patch_test_002.png :class: sphx-glr-multi-img * .. image-sg:: /auto_examples/plates/images/sphx_glr_plot_plate08_triangle_patch_test_003.png :alt: Deformed System (magnification=25.00) :srcset: /auto_examples/plates/images/sphx_glr_plot_plate08_triangle_patch_test_003.png :class: sphx-glr-multi-img * .. image-sg:: /auto_examples/plates/images/sphx_glr_plot_plate08_triangle_patch_test_004.png :alt: Contours of '$\sigma_{xx}$' :srcset: /auto_examples/plates/images/sphx_glr_plot_plate08_triangle_patch_test_004.png :class: sphx-glr-multi-img * .. image-sg:: /auto_examples/plates/images/sphx_glr_plot_plate08_triangle_patch_test_005.png :alt: Contours of '$\sigma_{yy}$' :srcset: /auto_examples/plates/images/sphx_glr_plot_plate08_triangle_patch_test_005.png :class: sphx-glr-multi-img * .. image-sg:: /auto_examples/plates/images/sphx_glr_plot_plate08_triangle_patch_test_006.png :alt: Contours of '$\sigma_{xy}$' :srcset: /auto_examples/plates/images/sphx_glr_plot_plate08_triangle_patch_test_006.png :class: sphx-glr-multi-img .. rst-class:: sphx-glr-script-out .. code-block:: none System Analysis Report ======================= Nodes: --------------------- Node_655: x: [0.000 0.000] fix: ['ux', 'uy'] u: None Node_656: x: [20.000 0.000] u: None Node_657: x: [50.000 0.000] u: None Node_658: x: [70.000 0.000] u: None Node_659: x: [100.000 0.000] fix: ['uy'] u: None Node_660: x: [0.000 16.000] u: None Node_661: x: [15.000 24.000] u: None Node_662: x: [50.000 16.000] u: None Node_663: x: [80.000 24.000] u: None Node_664: x: [100.000 16.000] u: None Node_665: x: [0.000 48.000] u: None Node_666: x: [20.000 40.000] u: None Node_667: x: [70.000 56.000] u: None Node_668: x: [90.000 48.000] u: None Node_669: x: [100.000 64.000] u: None Node_670: x: [0.000 80.000] u: None Node_671: x: [30.000 80.000] u: None Node_672: x: [55.000 80.000] u: None Node_673: x: [80.000 80.000] u: None Node_674: x: [100.000 80.000] u: None Elements: --------------------- Triangle_849: nodes ( Node_655 Node_656 Node_660 ) material: PlaneStress strain: xx=0.000e+00 yy=0.000e+00 xy=0.000e+00 zz=-0.000e+00 stress: xx=0.000e+00 yy=0.000e+00 xy=0.000e+00 zz=0.000e+00 Triangle_850: nodes ( Node_656 Node_657 Node_661 ) material: PlaneStress strain: xx=0.000e+00 yy=0.000e+00 xy=0.000e+00 zz=-0.000e+00 stress: xx=0.000e+00 yy=0.000e+00 xy=0.000e+00 zz=0.000e+00 Triangle_851: nodes ( Node_657 Node_658 Node_662 ) material: PlaneStress strain: xx=0.000e+00 yy=0.000e+00 xy=0.000e+00 zz=-0.000e+00 stress: xx=0.000e+00 yy=0.000e+00 xy=0.000e+00 zz=0.000e+00 Triangle_852: nodes ( Node_658 Node_659 Node_663 ) material: PlaneStress strain: xx=0.000e+00 yy=0.000e+00 xy=0.000e+00 zz=-0.000e+00 stress: xx=0.000e+00 yy=0.000e+00 xy=0.000e+00 zz=0.000e+00 Triangle_853: nodes ( Node_661 Node_660 Node_656 ) material: PlaneStress strain: xx=0.000e+00 yy=0.000e+00 xy=0.000e+00 zz=-0.000e+00 stress: xx=0.000e+00 yy=0.000e+00 xy=0.000e+00 zz=0.000e+00 Triangle_854: nodes ( Node_662 Node_661 Node_657 ) material: PlaneStress strain: xx=0.000e+00 yy=0.000e+00 xy=0.000e+00 zz=-0.000e+00 stress: xx=0.000e+00 yy=0.000e+00 xy=0.000e+00 zz=0.000e+00 Triangle_855: nodes ( Node_663 Node_662 Node_658 ) material: PlaneStress strain: xx=0.000e+00 yy=0.000e+00 xy=0.000e+00 zz=-0.000e+00 stress: xx=0.000e+00 yy=0.000e+00 xy=0.000e+00 zz=0.000e+00 Triangle_856: nodes ( Node_664 Node_663 Node_659 ) material: PlaneStress strain: xx=0.000e+00 yy=0.000e+00 xy=0.000e+00 zz=-0.000e+00 stress: xx=0.000e+00 yy=0.000e+00 xy=0.000e+00 zz=0.000e+00 Triangle_857: nodes ( Node_660 Node_661 Node_665 ) material: PlaneStress strain: xx=0.000e+00 yy=0.000e+00 xy=0.000e+00 zz=-0.000e+00 stress: xx=0.000e+00 yy=0.000e+00 xy=0.000e+00 zz=0.000e+00 Triangle_858: nodes ( Node_661 Node_662 Node_666 ) material: PlaneStress strain: xx=0.000e+00 yy=0.000e+00 xy=0.000e+00 zz=-0.000e+00 stress: xx=0.000e+00 yy=0.000e+00 xy=0.000e+00 zz=0.000e+00 Triangle_859: nodes ( Node_662 Node_663 Node_667 ) material: PlaneStress strain: xx=0.000e+00 yy=0.000e+00 xy=0.000e+00 zz=-0.000e+00 stress: xx=0.000e+00 yy=0.000e+00 xy=0.000e+00 zz=0.000e+00 Triangle_860: nodes ( Node_663 Node_664 Node_668 ) material: PlaneStress strain: xx=0.000e+00 yy=0.000e+00 xy=0.000e+00 zz=-0.000e+00 stress: xx=0.000e+00 yy=0.000e+00 xy=0.000e+00 zz=0.000e+00 Triangle_861: nodes ( Node_666 Node_665 Node_661 ) material: PlaneStress strain: xx=0.000e+00 yy=0.000e+00 xy=0.000e+00 zz=-0.000e+00 stress: xx=0.000e+00 yy=0.000e+00 xy=0.000e+00 zz=0.000e+00 Triangle_862: nodes ( Node_667 Node_666 Node_662 ) material: PlaneStress strain: xx=0.000e+00 yy=0.000e+00 xy=0.000e+00 zz=-0.000e+00 stress: xx=0.000e+00 yy=0.000e+00 xy=0.000e+00 zz=0.000e+00 Triangle_863: nodes ( Node_668 Node_667 Node_663 ) material: PlaneStress strain: xx=0.000e+00 yy=0.000e+00 xy=0.000e+00 zz=-0.000e+00 stress: xx=0.000e+00 yy=0.000e+00 xy=0.000e+00 zz=0.000e+00 Triangle_864: nodes ( Node_669 Node_668 Node_664 ) material: PlaneStress strain: xx=0.000e+00 yy=0.000e+00 xy=0.000e+00 zz=-0.000e+00 stress: xx=0.000e+00 yy=0.000e+00 xy=0.000e+00 zz=0.000e+00 Triangle_865: nodes ( Node_665 Node_666 Node_670 ) material: PlaneStress strain: xx=0.000e+00 yy=0.000e+00 xy=0.000e+00 zz=-0.000e+00 stress: xx=0.000e+00 yy=0.000e+00 xy=0.000e+00 zz=0.000e+00 Triangle_866: nodes ( Node_666 Node_667 Node_671 ) material: PlaneStress strain: xx=0.000e+00 yy=0.000e+00 xy=0.000e+00 zz=-0.000e+00 stress: xx=0.000e+00 yy=0.000e+00 xy=0.000e+00 zz=0.000e+00 Triangle_867: nodes ( Node_667 Node_668 Node_672 ) material: PlaneStress strain: xx=0.000e+00 yy=0.000e+00 xy=0.000e+00 zz=-0.000e+00 stress: xx=0.000e+00 yy=0.000e+00 xy=0.000e+00 zz=0.000e+00 Triangle_868: nodes ( Node_668 Node_669 Node_673 ) material: PlaneStress strain: xx=0.000e+00 yy=0.000e+00 xy=0.000e+00 zz=-0.000e+00 stress: xx=0.000e+00 yy=0.000e+00 xy=0.000e+00 zz=0.000e+00 Triangle_869: nodes ( Node_671 Node_670 Node_666 ) material: PlaneStress strain: xx=0.000e+00 yy=0.000e+00 xy=0.000e+00 zz=-0.000e+00 stress: xx=0.000e+00 yy=0.000e+00 xy=0.000e+00 zz=0.000e+00 Triangle_870: nodes ( Node_672 Node_671 Node_667 ) material: PlaneStress strain: xx=0.000e+00 yy=0.000e+00 xy=0.000e+00 zz=-0.000e+00 stress: xx=0.000e+00 yy=0.000e+00 xy=0.000e+00 zz=0.000e+00 Triangle_871: nodes ( Node_673 Node_672 Node_668 ) material: PlaneStress strain: xx=0.000e+00 yy=0.000e+00 xy=0.000e+00 zz=-0.000e+00 stress: xx=0.000e+00 yy=0.000e+00 xy=0.000e+00 zz=0.000e+00 Triangle_872: nodes ( Node_674 Node_673 Node_669 ) material: PlaneStress strain: xx=0.000e+00 yy=0.000e+00 xy=0.000e+00 zz=-0.000e+00 stress: xx=0.000e+00 yy=0.000e+00 xy=0.000e+00 zz=0.000e+00 + System Analysis Report ======================= Nodes: --------------------- Node_655: x: [0.000 0.000] fix: ['ux', 'uy'] u: [0.000 0.000] Node_656: x: [20.000 0.000] u: [0.100 0.000] Node_657: x: [50.000 0.000] u: [0.250 0.000] Node_658: x: [70.000 0.000] u: [0.350 -0.000] Node_659: x: [100.000 0.000] fix: ['uy'] u: [0.500 0.000] Node_660: x: [0.000 16.000] u: [-0.000 -0.020] Node_661: x: [15.000 24.000] u: [0.075 -0.030] Node_662: x: [50.000 16.000] u: [0.250 -0.020] Node_663: x: [80.000 24.000] u: [0.400 -0.030] Node_664: x: [100.000 16.000] u: [0.500 -0.020] Node_665: x: [0.000 48.000] u: [0.000 -0.060] Node_666: x: [20.000 40.000] u: [0.100 -0.050] Node_667: x: [70.000 56.000] u: [0.350 -0.070] Node_668: x: [90.000 48.000] u: [0.450 -0.060] Node_669: x: [100.000 64.000] u: [0.500 -0.080] Node_670: x: [0.000 80.000] u: [-0.000 -0.100] Node_671: x: [30.000 80.000] u: [0.150 -0.100] Node_672: x: [55.000 80.000] u: [0.275 -0.100] Node_673: x: [80.000 80.000] u: [0.400 -0.100] Node_674: x: [100.000 80.000] u: [0.500 -0.100] Elements: --------------------- Triangle_849: nodes ( Node_655 Node_656 Node_660 ) material: PlaneStress strain: xx=5.000e-03 yy=-1.250e-03 xy=-4.559e-17 zz=-9.375e-04 stress: xx=1.000e+02 yy=1.187e-12 xy=-3.647e-13 zz=0.000e+00 Triangle_850: nodes ( Node_656 Node_657 Node_661 ) material: PlaneStress strain: xx=5.000e-03 yy=-1.250e-03 xy=-7.546e-17 zz=-9.375e-04 stress: xx=1.000e+02 yy=0.000e+00 xy=-6.036e-13 zz=0.000e+00 Triangle_851: nodes ( Node_657 Node_658 Node_662 ) material: PlaneStress strain: xx=5.000e-03 yy=-1.250e-03 xy=-2.756e-16 zz=-9.375e-04 stress: xx=1.000e+02 yy=-1.183e-12 xy=-2.205e-12 zz=0.000e+00 Triangle_852: nodes ( Node_658 Node_659 Node_663 ) material: PlaneStress strain: xx=5.000e-03 yy=-1.250e-03 xy=-1.035e-16 zz=-9.375e-04 stress: xx=1.000e+02 yy=-1.183e-12 xy=-8.279e-13 zz=0.000e+00 Triangle_853: nodes ( Node_661 Node_660 Node_656 ) material: PlaneStress strain: xx=5.000e-03 yy=-1.250e-03 xy=2.776e-17 zz=-9.375e-04 stress: xx=1.000e+02 yy=3.553e-12 xy=2.220e-13 zz=0.000e+00 Triangle_854: nodes ( Node_662 Node_661 Node_657 ) material: PlaneStress strain: xx=5.000e-03 yy=-1.250e-03 xy=-5.551e-17 zz=-9.375e-04 stress: xx=1.000e+02 yy=-1.183e-12 xy=-4.441e-13 zz=0.000e+00 Triangle_855: nodes ( Node_663 Node_662 Node_658 ) material: PlaneStress strain: xx=5.000e-03 yy=-1.250e-03 xy=1.665e-16 zz=-9.375e-04 stress: xx=1.000e+02 yy=-1.183e-12 xy=1.332e-12 zz=0.000e+00 Triangle_856: nodes ( Node_664 Node_663 Node_659 ) material: PlaneStress strain: xx=5.000e-03 yy=-1.250e-03 xy=-5.551e-17 zz=-9.375e-04 stress: xx=1.000e+02 yy=1.183e-12 xy=-4.441e-13 zz=0.000e+00 Triangle_857: nodes ( Node_660 Node_661 Node_665 ) material: PlaneStress strain: xx=5.000e-03 yy=-1.250e-03 xy=3.163e-16 zz=-9.375e-04 stress: xx=1.000e+02 yy=-1.183e-12 xy=2.530e-12 zz=0.000e+00 Triangle_858: nodes ( Node_661 Node_662 Node_666 ) material: PlaneStress strain: xx=5.000e-03 yy=-1.250e-03 xy=5.551e-17 zz=-9.375e-04 stress: xx=1.000e+02 yy=-7.105e-12 xy=4.441e-13 zz=0.000e+00 Triangle_859: nodes ( Node_662 Node_663 Node_667 ) material: PlaneStress strain: xx=5.000e-03 yy=-1.250e-03 xy=1.110e-16 zz=-9.375e-04 stress: xx=1.000e+02 yy=4.736e-12 xy=8.882e-13 zz=0.000e+00 Triangle_860: nodes ( Node_663 Node_664 Node_668 ) material: PlaneStress strain: xx=5.000e-03 yy=-1.250e-03 xy=3.886e-16 zz=-9.375e-04 stress: xx=1.000e+02 yy=2.370e-12 xy=3.109e-12 zz=0.000e+00 Triangle_861: nodes ( Node_666 Node_665 Node_661 ) material: PlaneStress strain: xx=5.000e-03 yy=-1.250e-03 xy=-1.665e-16 zz=-9.375e-04 stress: xx=1.000e+02 yy=-8.288e-12 xy=-1.332e-12 zz=0.000e+00 Triangle_862: nodes ( Node_667 Node_666 Node_662 ) material: PlaneStress strain: xx=5.000e-03 yy=-1.250e-03 xy=5.551e-17 zz=-9.375e-04 stress: xx=1.000e+02 yy=0.000e+00 xy=4.441e-13 zz=0.000e+00 Triangle_863: nodes ( Node_668 Node_667 Node_663 ) material: PlaneStress strain: xx=5.000e-03 yy=-1.250e-03 xy=2.776e-16 zz=-9.375e-04 stress: xx=1.000e+02 yy=2.370e-12 xy=2.220e-12 zz=0.000e+00 Triangle_864: nodes ( Node_669 Node_668 Node_664 ) material: PlaneStress strain: xx=5.000e-03 yy=-1.250e-03 xy=1.395e-16 zz=-9.375e-04 stress: xx=1.000e+02 yy=-3.553e-12 xy=1.116e-12 zz=0.000e+00 Triangle_865: nodes ( Node_665 Node_666 Node_670 ) material: PlaneStress strain: xx=5.000e-03 yy=-1.250e-03 xy=-1.008e-16 zz=-9.375e-04 stress: xx=1.000e+02 yy=5.922e-12 xy=-8.065e-13 zz=0.000e+00 Triangle_866: nodes ( Node_666 Node_667 Node_671 ) material: PlaneStress strain: xx=5.000e-03 yy=-1.250e-03 xy=5.551e-17 zz=-9.375e-04 stress: xx=1.000e+02 yy=4.736e-12 xy=4.441e-13 zz=0.000e+00 Triangle_867: nodes ( Node_667 Node_668 Node_672 ) material: PlaneStress strain: xx=5.000e-03 yy=-1.250e-03 xy=1.110e-16 zz=-9.375e-04 stress: xx=1.000e+02 yy=-1.183e-12 xy=8.882e-13 zz=0.000e+00 Triangle_868: nodes ( Node_668 Node_669 Node_673 ) material: PlaneStress strain: xx=5.000e-03 yy=-1.250e-03 xy=-4.718e-16 zz=-9.375e-04 stress: xx=1.000e+02 yy=-3.553e-12 xy=-3.775e-12 zz=0.000e+00 Triangle_869: nodes ( Node_671 Node_670 Node_666 ) material: PlaneStress strain: xx=5.000e-03 yy=-1.250e-03 xy=-2.952e-17 zz=-9.375e-04 stress: xx=1.000e+02 yy=4.736e-12 xy=-2.362e-13 zz=0.000e+00 Triangle_870: nodes ( Node_672 Node_671 Node_667 ) material: PlaneStress strain: xx=5.000e-03 yy=-1.250e-03 xy=-1.268e-16 zz=-9.375e-04 stress: xx=1.000e+02 yy=-1.539e-11 xy=-1.014e-12 zz=0.000e+00 Triangle_871: nodes ( Node_673 Node_672 Node_668 ) material: PlaneStress strain: xx=5.000e-03 yy=-1.250e-03 xy=-4.441e-16 zz=-9.375e-04 stress: xx=1.000e+02 yy=-4.739e-12 xy=-3.553e-12 zz=0.000e+00 Triangle_872: nodes ( Node_674 Node_673 Node_669 ) material: PlaneStress strain: xx=5.000e-03 yy=-1.250e-03 xy=-0.000e+00 zz=-9.375e-04 stress: xx=1.000e+02 yy=-1.187e-12 xy=0.000e+00 zz=0.000e+00 /Users/pmackenz/Development/Educational/FEM.edu/venv/lib/python3.13/site-packages/matplotlib/quiver.py:678: RuntimeWarning: Mean of empty slice. amean = a.mean() /Users/pmackenz/Development/Educational/FEM.edu/venv/lib/python3.13/site-packages/numpy/_core/_methods.py:145: RuntimeWarning: invalid value encountered in scalar divide ret = ret.dtype.type(ret / rcount) .. rst-class:: sphx-glr-timing **Total running time of the script:** (0 minutes 0.161 seconds) .. _sphx_glr_download_auto_examples_plates_plot_plate08_triangle_patch_test.py: .. only:: html .. container:: sphx-glr-footer sphx-glr-footer-example .. container:: sphx-glr-download sphx-glr-download-jupyter :download:`Download Jupyter notebook: plot_plate08_triangle_patch_test.ipynb ` .. container:: sphx-glr-download sphx-glr-download-python :download:`Download Python source code: plot_plate08_triangle_patch_test.py ` .. container:: sphx-glr-download sphx-glr-download-zip :download:`Download zipped: plot_plate08_triangle_patch_test.zip ` .. only:: html .. rst-class:: sphx-glr-signature `Gallery generated by Sphinx-Gallery `_