我有一个二维数据和四边形的二维网格描述一个细分为补丁的域.数据在每个网格节点处定义.数据的不连续性存在于补丁边界处,即数据在同一位置被多次定义.
如何使用Python在节点之间使用线性插值绘制此数据,并在每个补丁面上正确表示不连续值?
下面是三个示例元素或补丁,每个元素或补丁各有六个节点值.

节点位置和值数据可以存储在[Kx3x2]数组中,其中K是元素的数量.例如,
x = np.array( [
[ [0.0, 1.0], [0.0, 1.0], [0.0, 1.0] ], #element 0
[ [1.0, 2.0], [1.0, 2.0], [1.0, 2.0] ], #element 1
[ [2.0, 3.0], [2.0, 3.0], [2.0, 3.0] ], #element 2
] )
y = np.array( [
[ [0.0, 0.0], [0.5, 0.5], [1.0, 1.0] ], #element 0
[ [0.0, 1.0], [0.5, 1.5], [1.0, 2.0] ], #element 1
[ [1.0, 1.0], [1.5, 1.5], [2.0, 2.0] ], #element 2
] )
z = np.array( [
[ [0.0, 0.5], [0.0, 0.8], [0.0, 1.0] ], #element 0
[ [0.3, 1.0], [0.6, 1.2], [0.8, 1.3] ], #element 1
[ [1.2, 1.5], [1.3, 1.4], [1.5, 1.7] ], #element 2
] )
Run Code Online (Sandbox Code Playgroud)
我考虑过了pyplot.imshow().这不能同时考虑整个域,并且仍然代表多值不连续节点.可能需要imshow()为每个补丁单独调用.但是,我如何在同一轴上绘制每个补丁图像? imshow()对于非矩形补丁也是有问题的,这是我的一般情况.
我考虑过pyplot.pcolormesh(),但它似乎只适用于以细胞为中心的数据.
一个选项通过对所有元素进行三角测量,然后使用tripcolor()我现在发现的matplotlib 函数进行绘图.这里和这里有两个有用的演示.
我的全局域的自动三角测量可能会有问题,但单个四边形的Delaunay三角剖分非常有效:

我通过附加每个元素的三角测量来创建全局三角测量.这意味着共享节点实际上在位置数组和值数组中重复.这允许元素面处的不连续数据.

可以使用tripcolor()函数,为每个节点提供节点位置和值来实现根据需要绘制线性插值和不连续性.

我有点担心轮廓绘图可能如何工作,因为元素面不再逻辑连接. tricontour()仍然按预期工作.(此处显示三角形覆盖)

使用以下代码转载:
import numpy as np
import matplotlib.pyplot as plt
import matplotlib.tri as tri
x = np.array( [
[ [0.0, 1.0], [0.0, 1.0], [0.0, 1.0] ], #element 0
[ [1.0, 2.0], [1.0, 2.0], [1.0, 2.0] ], #element 1
[ [2.0, 3.0], [2.0, 3.0], [2.0, 3.0] ], #element 2
] )
y = np.array( [
[ [0.0, 0.0], [0.5, 0.5], [1.0, 1.0] ], #element 0
[ [0.0, 1.0], [0.5, 1.5], [1.0, 2.0] ], #element 1
[ [1.0, 1.0], [1.5, 1.5], [2.0, 2.0] ], #element 2
] )
z = np.array( [
[ [0.0, 0.5], [0.0, 0.8], [0.0, 1.0] ], #element 0
[ [0.3, 1.0], [0.6, 1.2], [0.8, 1.3] ], #element 1
[ [1.2, 1.5], [1.3, 1.4], [1.5, 1.7] ], #element 2
] )
global_num_pts = z.size
global_x = np.zeros( global_num_pts )
global_y = np.zeros( global_num_pts )
global_z = np.zeros( global_num_pts )
global_triang_list = list()
offset = 0;
num_triangles = 0;
#process triangulation element-by-element
for k in range(z.shape[0]):
points_x = x[k,...].flatten()
points_y = y[k,...].flatten()
z_element = z[k,...].flatten()
num_points_this_element = points_x.size
#auto-generate Delauny triangulation for the element, which should be flawless due to quadrilateral element shape
triang = tri.Triangulation(points_x, points_y)
global_triang_list.append( triang.triangles + offset ) #offseting triangle indices by start index of this element
#store results for this element in global triangulation arrays
global_x[offset:(offset+num_points_this_element)] = points_x
global_y[offset:(offset+num_points_this_element)] = points_y
global_z[offset:(offset+num_points_this_element)] = z_element
num_triangles += triang.triangles.shape[0]
offset += num_points_this_element
#go back and turn all of the triangle indices into one global triangle array
offset = 0
global_triang = np.zeros( (num_triangles, 3) )
for t in global_triang_list:
global_triang[ offset:(offset+t.shape[0] )] = t
offset += t.shape[0]
plt.figure()
plt.gca().set_aspect('equal')
plt.tripcolor(global_x, global_y, global_triang, global_z, shading='gouraud' )
#plt.tricontour(global_x, global_y, global_triang, global_z )
#plt.triplot(global_x, global_y, global_triang, 'go-') #plot just the triangle mesh
plt.xlim((-0.25, 3.25))
plt.ylim((-0.25, 2.25))
plt.show()
Run Code Online (Sandbox Code Playgroud)