Wrong plot with plot_trisurf

General
1

I am trying to make a surface plot of a solution on posed on a domain which is a nonconvex polygon, specifically, the square [-1,1]x[-1,1] minus the square [0,1]x[0,1]. I am using this code (which seems to be too long!)

#####
import matplotlib.pyplot as plt
from mpl_toolkits.mplot3d import Axes3D
#####
#####

# graphics
  coordinates = mesh2d.coordinates()
  XX = [0] * (3*nel[k])
  YY = [0] * (3*nel[k])
  ZZ = [0] * (3*nel[k])
  dofmap = V.dofmap()
  nodal_values =  u.vector()
  dof2coord= coordinates[dof_to_vertex_map(V)]
  for cell in cells(mesh2d):
    dof = dofmap.cell_dofs(cell.index())
    x0 = dof2coord[dof[0],0]
    x1 = dof2coord[dof[1],0]
    x2 = dof2coord[dof[2],0]
    y0 = dof2coord[dof[0],1]
    y1 = dof2coord[dof[1],1]
    y2 = dof2coord[dof[2],1]
    #check orientation
    if (x1-x0)*(y2-y0)-(x2-x0)*(y1-y0)>0:
      i0 = 0
      i1 = 1
      i2 = 2
    else:
      i0 = 0
      i1 = 2
      i2 = 1
    XX[3*cell.index()] = dof2coord[dof[i0],0]
    XX[3*cell.index()+1] = dof2coord[dof[i1],0]
    XX[3*cell.index()+2] = dof2coord[dof[i2],0]
    YY[3*cell.index()] = dof2coord[dof[i0],1]
    YY[3*cell.index()+1] = dof2coord[dof[i1],1]
    YY[3*cell.index()+2] = dof2coord[dof[i2],1]
    ZZ[3*cell.index()] = nodal_values[dof[i0]]
    ZZ[3*cell.index()+1] = nodal_values[dof[i1]]
    ZZ[3*cell.index()+2] = nodal_values[dof[i2]]  

  fig = plt.figure(figsize =(16, 9))
  ax = fig.add_subplot(projection='3d')
  # Creating color map
  my_cmap = plt.get_cmap('summer')
  trisurf = ax.plot_trisurf(XX, YY, ZZ, linewidth=0.2, antialiased=False, cmap=my_cmap)
  fig.colorbar(trisurf, ax = ax, shrink = 0.5, aspect = 5)
  plt.show()

It results that it draws some triangles which do not correspond to the plot, I mean, close to the nonconvex vertex there appear wrong triangles, like the ones showed in this figure:

Captura
Captura779×638 158 KB

I wanted to know if I am using in a wrong way the command plot_trisurf. Thanks!