如何使用matplotlib绘制二维有限元结果？

Question

我正在开发2D平面有限元工具。功能之一是能够可视化特定对象上的压力。

此工具使用以下数据创建四边形网格：

• 节点：numpy数组[[x1 y1], [x2 y2], etc]->网格中每个节点的x和y坐标

• elements：numpy数组[[1 2 3 4], [2 3 5 6]]->数组的每一行对应于网格的一个特定元素的4个点。

我能够实现绘制网格的方法：

import matplotlib.pyplot as plt
import matplotlib.collections
import matplotlib.cm as cm

import numpy as np


def showMeshPlot(nodes, elements):

    y = nodes[:,0]
    z = nodes[:,1]

    #https://stackoverflow.com/questions/49640311/matplotlib-unstructered-quadrilaterals-instead-of-triangles
    def quatplot(y,z, quatrangles, ax=None, **kwargs):

        if not ax: ax=plt.gca()
        yz = np.c_[y,z]
        verts= yz[quatrangles]
        pc = matplotlib.collections.PolyCollection(verts, **kwargs)
        ax.add_collection(pc)
        ax.autoscale()

    plt.figure()
    plt.gca().set_aspect('equal')

    quatplot(y,z, np.asarray(elements), ax=None, color="crimson", facecolor="None")
    if nodes:            
        plt.plot(y,z, marker="o", ls="", color="crimson")

    plt.title('This is the plot for: quad')
    plt.xlabel('Y Axis')
    plt.ylabel('Z Axis')


    plt.show()

nodes = np.array([[0,0], [0,0.5],[0,1],[0.5,0], [0.5,0.5], [0.5,1], [1,0], 
                  [1,0.5],[1,1]])
elements = np.array([[0,3,4,1],[1,4,5,2],[3,6,7,4],[4,7,8,5]])
stresses = np.array([1,2,3,4])

showMeshPlot(nodes, elements)

哪个会产生这样的情节：

现在，我有一个一维数组，对象上的应力与元素数组的长度相同。

我的问题是如何使用matplotlib可视化那些应力（使用标尺）？我调查了pcolormesh，但是我不明白它如何处理我的数据。这是我要达到的目标的一个示例（对robbievanleeuwen的积分）：

注意：我无法复制上面的示例，因为他使用了三角形网格而不是四边形。

提前感谢！

Answer 1

PolyCollection是ScalarMappable。它可以具有值数组，颜色图和规范化集。在这里，您可以将stresses数组提供给PolyCollection并选择一些要使用的颜色图。其余的重新调整了该功能，以便可以将其他数据作为输入并创建一个颜色条。

import matplotlib.pyplot as plt
import matplotlib.collections
import numpy as np


def showMeshPlot(nodes, elements, values):

    y = nodes[:,0]
    z = nodes[:,1]

    def quatplot(y,z, quatrangles, values, ax=None, **kwargs):

        if not ax: ax=plt.gca()
        yz = np.c_[y,z]
        verts= yz[quatrangles]
        pc = matplotlib.collections.PolyCollection(verts, **kwargs)
        pc.set_array(values)
        ax.add_collection(pc)
        ax.autoscale()
        return pc

    fig, ax = plt.subplots()
    ax.set_aspect('equal')

    pc = quatplot(y,z, np.asarray(elements), values, ax=ax, 
             edgecolor="crimson", cmap="rainbow")
    fig.colorbar(pc, ax=ax)        
    ax.plot(y,z, marker="o", ls="", color="crimson")

    ax.set(title='This is the plot for: quad', xlabel='Y Axis', ylabel='Z Axis')

    plt.show()

nodes = np.array([[0,0], [0,0.5],[0,1],[0.5,0], [0.5,0.5], [0.5,1], [1,0], 
                  [1,0.5],[1,1]])
elements = np.array([[0,3,4,1],[1,4,5,2],[3,6,7,4],[4,7,8,5]])
stresses = np.array([1,2,3,4])

showMeshPlot(nodes, elements, stresses)

Answer 2

我认为您最好的选择是使用tricontour。您已经有三角剖分了，对吧？

它创建了这样的图：

（[from here）

https://matplotlib.org/api/_as_gen/matplotlib.axes.Axes.tricontour.html

Answer 3

您可以使用matplotlib.pyplot.tricontourf(x, y, triangles, component, ticks, **kwargs)：

[x和y是一维列表或包含网格节点坐标的numpy.arrays。

triangles是一个number_of_triangles除以3（每个三角形具有3个节点）的2D列表或numpy.array。例如，如果此数组的第一行是[0，1，2]，则第一个三角形将连接节点0、1和2。请记住，在python数组中，数组以0条目开头。

[component是一维列表或numpy.array，其中包含要绘制的组件的节点值，例如节点上的von Mises应力。

[ticks定义所需的轮廓，可以定义ticks，例如：

ticks = numpy.linspace(min(component), max(component), 13) # this will give you 12 colors and 13 ticks on the colorbar

一旦调用该函数，您将获得如下内容：

请注意，某些元素是四边形，为此，每个四边形必须拆分为2个三角形。因此，如果您有一个包含4个四边形的网格，则可以将它们分成8个三角形，如下所示：

def QuadsToTris(quads):
    triquads = [[None for c in range(3)] for l in range(2*len(quads))]
    for i in range(0, len(quads), 1):
        j = 2*i
        n1 = quads[i][0]
        n2 = quads[i][1]
        n3 = quads[i][2]
        n4 = quads[i][3]
        triquads[j][0] = n1; triquads[j+1][0] = n3
        triquads[j][1] = n2; triquads[j+1][1] = n4
        triquads[j][2] = n3; triquads[j+1][2] = n1
    return triquads #this list contains the triangles made by splitting the quadrangles

然后只绘制带有新三角形的轮廓。

如果您实际上想看到四边形网格（黑线），则在绘制轮廓后，可以使用matplotlib.pyplot.plot（x，y）]绘制每个单独的四边形：>]

def PlotQuadMesh(quads, vx, vy): #vx and vy are the previous x and y on tricontourf
    for i in range(0, len(quads), 1):
        x = [None]*5; y = [None]*5
        n1 = quads[i][0]
        n2 = quads[i][1]
        n3 = quads[i][2]
        n4 = quads[i][3]
        x[0] = vx[n1]; y[0] = vy[n1]
        x[1] = vx[n2]; y[1] = vy[n2]
        x[2] = vx[n3]; y[2] = vy[n3]
        x[3] = vx[n4]; y[3] = vy[n4]
        x[4] = vx[n1]; y[4] = vy[n1] #don't forget to close the element
        plt.plot(x, y, color = 'black', linestyle = '-', linewidth = '0.5')
这种绘制网格的方法缓慢且效率低下，请考虑使用VTK。