diff --git a/python/core/auto_generated/mesh/qgsmeshdataprovider.sip.in b/python/core/auto_generated/mesh/qgsmeshdataprovider.sip.in
index 34cecf2a30a..5eafaf5f8cd 100644
--- a/python/core/auto_generated/mesh/qgsmeshdataprovider.sip.in
+++ b/python/core/auto_generated/mesh/qgsmeshdataprovider.sip.in
@@ -84,6 +84,8 @@ Compare two faces, return ``True`` if they are equivalent : same indexes and sam
 
 };
 
+
+
 class QgsMeshDataSourceInterface /Abstract/
 {
 %Docstring(signature="appended")
diff --git a/python/core/auto_generated/mesh/qgsmeshlayer.sip.in b/python/core/auto_generated/mesh/qgsmeshlayer.sip.in
index 1598ed3d485..ad3f120c534 100644
--- a/python/core/auto_generated/mesh/qgsmeshlayer.sip.in
+++ b/python/core/auto_generated/mesh/qgsmeshlayer.sip.in
@@ -917,7 +917,7 @@ Access to labeling configuration. May be ``None`` if labeling is not used.
 
 .. seealso:: :py:func:`labelsEnabled`
 
-.. versionadded:: 3.26
+.. versionadded:: 3.36
 %End
 
     void setLabeling( QgsAbstractMeshLayerLabeling *labeling /Transfer/ );
diff --git a/python/core/auto_generated/mesh/qgsmeshlayerlabeling.sip.in b/python/core/auto_generated/mesh/qgsmeshlayerlabeling.sip.in
index 4e5d51099b7..d9c0f7c945c 100644
--- a/python/core/auto_generated/mesh/qgsmeshlayerlabeling.sip.in
+++ b/python/core/auto_generated/mesh/qgsmeshlayerlabeling.sip.in
@@ -18,7 +18,7 @@ class QgsAbstractMeshLayerLabeling
 %Docstring(signature="appended")
 Abstract base class - its implementations define different approaches to the labeling of a mesh layer.
 
-.. versionadded:: 3.38
+.. versionadded:: 3.36
 %End
 
 %TypeHeaderCode
@@ -65,8 +65,6 @@ Set pal settings for a specific provider (takes ownership).
 
 :param settings: Pal layer settings
 :param providerId: The id of the provider
-
-.. versionadded:: 3.0
 %End
 
     virtual bool requiresAdvancedEffects() const = 0;
@@ -82,8 +80,6 @@ Multiply opacity by ``opacityFactor``.
 
 This method multiplies the opacity of the labeling elements (text, shadow, buffer etc.)
 by ``opacity`` effectively changing the opacity of the whole labeling elements.
-
-.. versionadded:: 3.32
 %End
 
 
@@ -120,16 +116,18 @@ class QgsMeshLayerSimpleLabeling : QgsAbstractMeshLayerLabeling
 %Docstring(signature="appended")
 Basic implementation of the labeling interface for mesh layer.
 
-.. versionadded:: 3.38
+.. versionadded:: 3.36
 %End
 
 %TypeHeaderCode
 #include "qgsmeshlayerlabeling.h"
 %End
   public:
-    explicit QgsMeshLayerSimpleLabeling( const QgsPalLayerSettings &settings );
+
+    explicit QgsMeshLayerSimpleLabeling( const QgsPalLayerSettings &settings, bool labelFaces = false );
 %Docstring
-Constructs simple labeling configuration with given initial settings
+Constructs simple labeling configuration with given initial ``settings``.
+Labels are placed on mesh vertices unless ``labelFaces`` is ``True``, when they are placed on mesh faces.
 %End
 
     virtual QString type() const;
diff --git a/resources/function_help/json/mesh_contour b/resources/function_help/json/mesh_contour
deleted file mode 100644
index 92501b7d915..00000000000
--- a/resources/function_help/json/mesh_contour
+++ /dev/null
@@ -1,24 +0,0 @@
-{
-  "name": "mesh_contour",
-  "type": "function",
-  "groups": ["Meshes"],
-  "description": "Returns the mesh scalar value at a given point for a dataset group used for contour rendering.",
-  "arguments": [{
-    "arg": "point",
-    "description": "point geometry (for multipart geometries having more than one part, a first part will be used)."
-  }, {
-    "arg": "timestamp",
-    "description": "timestamp (defaults to current time)."
-  }, {
-    "arg": "layer",
-    "description": "mesh layer (defaults to current layer)."
-  }],
-  "examples": [{
-    "expression": "mesh_contour(make_point(1,1))",
-    "returns": "2.5"
-  }, {
-    "expression": "mesh_contour(make_point(1,1), make_datetime(2020,5,4,13,45,30.5))",
-    "returns": "2.5"
-  }],
-  "tags": ["mesh", "contour", "point"]
-}
diff --git a/resources/function_help/json/mesh_data b/resources/function_help/json/mesh_data
deleted file mode 100644
index 8e5aaf8a4d2..00000000000
--- a/resources/function_help/json/mesh_data
+++ /dev/null
@@ -1,27 +0,0 @@
-{
-  "name": "mesh_data",
-  "type": "function",
-  "groups": ["Meshes"],
-  "description": "Returns the mesh value at a given point for a given dataset group",
-  "arguments": [{
-    "arg": "point",
-    "description": "point geometry (for multipart geometries having more than one part, a first part will be used)."
-  }, {
-    "arg": "dataset_group",
-    "description": "name of the mesh dataset group."
-  }, {
-    "arg": "timestamp",
-    "description": "timestamp (defaults to current time)."
-  }, {
-    "arg": "layer",
-    "description": "mesh layer (defaults to current layer)."
-  }],
-  "examples": [{
-    "expression": "mesh_data(make_point(1,1), 'Bed Elevation')",
-    "returns": "2.5"
-  }, {
-    "expression": "mesh_contour(make_point(1,1), 'Bed Elevation' make_datetime(2020,5,4,13,45,30.5))",
-    "returns": "2.5"
-  }],
-  "tags": ["mesh", "dataset", "point"]
-}
diff --git a/src/core/expression/qgsexpressioncontextutils.cpp b/src/core/expression/qgsexpressioncontextutils.cpp
index ba309c3421f..49e63b331f4 100644
--- a/src/core/expression/qgsexpressioncontextutils.cpp
+++ b/src/core/expression/qgsexpressioncontextutils.cpp
@@ -1043,6 +1043,17 @@ class CurrentVertexZValueExpressionFunction: public QgsScopedExpressionFunction
       if ( !context )
         return QVariant();
 
+      if ( context->hasVariable( QStringLiteral( "_mesh_vertex_index" ) ) && context->hasVariable( QStringLiteral( "_native_mesh" ) ) )
+      {
+        int vertexIndex = context->variable( QStringLiteral( "_mesh_vertex_index" ) ).toInt();
+        QgsMesh *nativeMesh = qvariant_cast<QgsMesh *>( context->variable( QStringLiteral( "_native_mesh" ) ) );
+        const QgsMeshVertex &vertex = nativeMesh->vertex( vertexIndex );
+        if ( !vertex.isEmpty() )
+          return vertex.z();
+        else
+          return QVariant();
+      }
+
       if ( !context->hasVariable( QStringLiteral( "_mesh_vertex_index" ) ) || !context->hasVariable( QStringLiteral( "_mesh_layer" ) ) )
         return QVariant();
 
@@ -1081,6 +1092,17 @@ class CurrentVertexXValueExpressionFunction: public QgsScopedExpressionFunction
       if ( !context )
         return QVariant();
 
+      if ( context->hasVariable( QStringLiteral( "_mesh_vertex_index" ) ) && context->hasVariable( QStringLiteral( "_native_mesh" ) ) )
+      {
+        int vertexIndex = context->variable( QStringLiteral( "_mesh_vertex_index" ) ).toInt();
+        QgsMesh *nativeMesh = qvariant_cast<QgsMesh *>( context->variable( QStringLiteral( "_native_mesh" ) ) );
+        const QgsMeshVertex &vertex = nativeMesh->vertex( vertexIndex );
+        if ( !vertex.isEmpty() )
+          return vertex.x();
+        else
+          return QVariant();
+      }
+
       if ( !context->hasVariable( QStringLiteral( "_mesh_vertex_index" ) ) || !context->hasVariable( QStringLiteral( "_mesh_layer" ) ) )
         return QVariant();
 
@@ -1119,6 +1141,17 @@ class CurrentVertexYValueExpressionFunction: public QgsScopedExpressionFunction
       if ( !context )
         return QVariant();
 
+      if ( context->hasVariable( QStringLiteral( "_mesh_vertex_index" ) ) && context->hasVariable( QStringLiteral( "_native_mesh" ) ) )
+      {
+        int vertexIndex = context->variable( QStringLiteral( "_mesh_vertex_index" ) ).toInt();
+        QgsMesh *nativeMesh = qvariant_cast<QgsMesh *>( context->variable( QStringLiteral( "_native_mesh" ) ) );
+        const QgsMeshVertex &vertex = nativeMesh->vertex( vertexIndex );
+        if ( !vertex.isEmpty() )
+          return vertex.y();
+        else
+          return QVariant();
+      }
+
       if ( !context->hasVariable( QStringLiteral( "_mesh_vertex_index" ) ) || !context->hasVariable( QStringLiteral( "_mesh_layer" ) ) )
         return QVariant();
 
@@ -1157,6 +1190,17 @@ class CurrentVertexExpressionFunction: public QgsScopedExpressionFunction
       if ( !context )
         return QVariant();
 
+      if ( context->hasVariable( QStringLiteral( "_mesh_vertex_index" ) ) && context->hasVariable( QStringLiteral( "_native_mesh" ) ) )
+      {
+        int vertexIndex = context->variable( QStringLiteral( "_mesh_vertex_index" ) ).toInt();
+        QgsMesh *nativeMesh = qvariant_cast<QgsMesh *>( context->variable( QStringLiteral( "_native_mesh" ) ) );
+        const QgsMeshVertex &vertex = nativeMesh->vertex( vertexIndex );
+        if ( !vertex.isEmpty() )
+          return QVariant::fromValue( QgsGeometry( new QgsPoint( vertex ) ) );
+        else
+          return QVariant();
+      }
+
       if ( !context->hasVariable( QStringLiteral( "_mesh_vertex_index" ) ) || !context->hasVariable( QStringLiteral( "_mesh_layer" ) ) )
         return QVariant();
 
@@ -1195,7 +1239,7 @@ class CurrentVertexIndexExpressionFunction: public QgsScopedExpressionFunction
       if ( !context )
         return QVariant();
 
-      if ( !context->hasVariable( QStringLiteral( "_mesh_vertex_index" ) ) || !context->hasVariable( QStringLiteral( "_mesh_layer" ) ) )
+      if ( !context->hasVariable( QStringLiteral( "_mesh_vertex_index" ) ) )
         return QVariant();
 
       return context->variable( QStringLiteral( "_mesh_vertex_index" ) );
@@ -1224,6 +1268,32 @@ class CurrentFaceAreaExpressionFunction: public QgsScopedExpressionFunction
       if ( !context )
         return QVariant();
 
+      if ( context->hasVariable( QStringLiteral( "_mesh_face_index" ) ) && context->hasVariable( QStringLiteral( "_native_mesh" ) ) )
+      {
+        const int faceIndex = context->variable( QStringLiteral( "_mesh_face_index" ) ).toInt();
+        QgsMesh *nativeMesh = qvariant_cast<QgsMesh *>( context->variable( QStringLiteral( "_native_mesh" ) ) );
+        const QgsMeshFace &face = nativeMesh->face( faceIndex );
+        if ( !face.isEmpty() )
+        {
+          QgsDistanceArea *calc = parent->geomCalculator();
+          QgsGeometry geom = QgsMeshUtils::toGeometry( nativeMesh->face( faceIndex ), nativeMesh->vertices );
+          if ( calc )
+          {
+            double area = calc->measureArea( geom );
+            area = calc->convertAreaMeasurement( area, parent->areaUnits() );
+            return QVariant( area );
+          }
+          else
+          {
+            return QVariant( geom.area() );
+          }
+        }
+        else
+        {
+          return QVariant();
+        }
+      }
+
       if ( !context->hasVariable( QStringLiteral( "_mesh_face_index" ) ) || !context->hasVariable( QStringLiteral( "_mesh_layer" ) ) )
         return QVariant();
 
@@ -1275,7 +1345,7 @@ class CurrentFaceIndexExpressionFunction: public QgsScopedExpressionFunction
       if ( !context )
         return QVariant();
 
-      if ( !context->hasVariable( QStringLiteral( "_mesh_face_index" ) ) || !context->hasVariable( QStringLiteral( "_mesh_layer" ) ) )
+      if ( !context->hasVariable( QStringLiteral( "_mesh_face_index" ) ) )
         return QVariant();
 
       return context->variable( QStringLiteral( "_mesh_face_index" ) ).toInt();
diff --git a/src/core/expression/qgsexpressionfunction.cpp b/src/core/expression/qgsexpressionfunction.cpp
index e150c67fe64..2f0fdcf3adb 100644
--- a/src/core/expression/qgsexpressionfunction.cpp
+++ b/src/core/expression/qgsexpressionfunction.cpp
@@ -64,8 +64,6 @@
 #include "qgsunittypes.h"
 #include "qgsspatialindex.h"
 #include "qgscolorrampimpl.h"
-#include "qgsmeshlayer.h"
-#include "qgsmeshdataset.h"
 
 #include <QMimeDatabase>
 #include <QProcessEnvironment>
@@ -1867,173 +1865,6 @@ static QVariant fcnRasterAttributes( const QVariantList &values, const QgsExpres
   }
 }
 
-static QVariant fcnMeshContour( const QVariantList &values, const QgsExpressionContext *context, QgsExpression *parent, const QgsExpressionNodeFunction * )
-{
-  const QgsGeometry geom = QgsExpressionUtils::getGeometry( values.at( 0 ), parent );
-
-  QDateTime datetime;
-  if ( values.size() < 2 || QgsVariantUtils::isNull( values.at( 1 ) ) )
-  {
-    datetime = QDateTime::currentDateTimeUtc();
-  }
-  else
-  {
-    datetime = QgsExpressionUtils::getDateTimeValue( values.at( 1 ), parent );
-  }
-
-  QVariant layer;
-  if ( values.size() < 3 || QgsVariantUtils::isNull( values.at( 2 ) ) )
-  {
-    layer = context->variable( QStringLiteral( "layer" ) );
-  }
-  else
-  {
-    layer = values.at( 2 );
-  }
-
-  bool foundLayer = false;
-  const QVariant res = QgsExpressionUtils::runMapLayerFunctionThreadSafe( layer, context, parent, [geom, parent]( QgsMapLayer * mapLayer ) -> QVariant
-  {
-    QgsMeshLayer *layer = qobject_cast< QgsMeshLayer * >( mapLayer );
-    if ( !layer )
-    {
-      return QVariant();
-    }
-
-    if ( geom.isNull() || geom.type() != Qgis::GeometryType::Point )
-    {
-      parent->setEvalErrorString( QObject::tr( "Function `mesh_contour` requires a valid point geometry." ) );
-      return QVariant();
-    }
-
-    QgsPointXY point = geom.asPoint();
-    if ( geom.isMultipart() )
-    {
-      QgsMultiPointXY multiPoint = geom.asMultiPoint();
-      if ( multiPoint.count() == 1 )
-      {
-        point = multiPoint[0];
-      }
-      else
-      {
-        return QVariant();
-      }
-    }
-
-    QgsMeshDatasetIndex index = layer->staticScalarDatasetIndex();
-    if ( !index.isValid() )
-    {
-      return QVariant();
-    }
-    const QgsMeshDatasetValue scalarValue = layer->datasetValue( index, point );
-    return scalarValue.scalar();
-  }, foundLayer );
-
-
-  if ( !foundLayer )
-  {
-    return QVariant();
-  }
-  else
-  {
-    return res;
-  }
-}
-
-static QVariant fcnMeshData( const QVariantList &values, const QgsExpressionContext *context, QgsExpression *parent, const QgsExpressionNodeFunction * )
-{
-  const QgsGeometry geom = QgsExpressionUtils::getGeometry( values.at( 0 ), parent );
-
-  QString datasetGroupName = values.at( 1 ).toString();
-  if ( datasetGroupName.isEmpty() )
-  {
-    parent->setEvalErrorString( QObject::tr( "Mesh dataset group name can not be empty." ) );
-    return QVariant();
-  }
-
-  QDateTime datetime;
-  if ( values.size() < 3 || QgsVariantUtils::isNull( values.at( 2 ) ) )
-  {
-    datetime = QDateTime::currentDateTimeUtc();
-  }
-  else
-  {
-    datetime = QgsExpressionUtils::getDateTimeValue( values.at( 2 ), parent );
-  }
-
-  QVariant layer;
-  if ( values.size() < 4 || QgsVariantUtils::isNull( values.at( 3 ) ) )
-  {
-    layer = context->variable( QStringLiteral( "layer" ) );
-  }
-  else
-  {
-    layer = values.at( 3 );
-  }
-
-  bool foundLayer = false;
-  const QVariant res = QgsExpressionUtils::runMapLayerFunctionThreadSafe( layer, context, parent, [geom, datasetGroupName, parent]( QgsMapLayer * mapLayer ) -> QVariant
-  {
-    QgsMeshLayer *layer = qobject_cast< QgsMeshLayer * >( mapLayer );
-    if ( !layer )
-    {
-      return QVariant();
-    }
-
-    if ( geom.isNull() || geom.type() != Qgis::GeometryType::Point )
-    {
-      parent->setEvalErrorString( QObject::tr( "Function `mesh_data` requires a valid point geometry." ) );
-      return QVariant();
-    }
-
-    QgsPointXY point = geom.asPoint();
-    if ( geom.isMultipart() )
-    {
-      QgsMultiPointXY multiPoint = geom.asMultiPoint();
-      if ( multiPoint.count() == 1 )
-      {
-        point = multiPoint[0];
-      }
-      else
-      {
-        return QVariant();
-      }
-    }
-
-    QgsMeshDatasetGroupTreeItem *root = layer->datasetGroupTreeRootItem();
-    QList<QgsMeshDatasetIndex> datasetIndexList;
-    const QList<int> allGroup = layer->enabledDatasetGroupsIndexes();
-    for ( int groupIndex : allGroup )
-    {
-      QgsMeshDatasetGroupTreeItem *group = root->childFromDatasetGroupIndex( groupIndex );
-      if ( group->name() == datasetGroupName )
-      {
-        datasetIndexList.append( QgsMeshDatasetIndex( groupIndex, 0 ) );
-        break;
-      }
-    }
-
-    if ( datasetIndexList.size() == 0 )
-    {
-      parent->setEvalErrorString( QObject::tr( "Dataset group '%1' not found." ).arg( datasetGroupName ) );
-      return QVariant();
-    }
-
-    const QgsMeshDatasetValue scalarValue = layer->datasetValue( datasetIndexList.at( 0 ), point );
-    return scalarValue.scalar();
-  }, foundLayer );
-
-
-  if ( !foundLayer )
-  {
-    return QVariant();
-  }
-  else
-  {
-    return res;
-  }
-}
-
 static QVariant fcnFeature( const QVariantList &, const QgsExpressionContext *context, QgsExpression *, const QgsExpressionNodeFunction * )
 {
   if ( !context )
@@ -9364,10 +9195,6 @@ const QList<QgsExpressionFunction *> &QgsExpression::Functions()
         << new QgsStaticExpressionFunction( QStringLiteral( "raster_value" ), QgsExpressionFunction::ParameterList() << QgsExpressionFunction::Parameter( QStringLiteral( "layer" ) ) << QgsExpressionFunction::Parameter( QStringLiteral( "band" ) ) << QgsExpressionFunction::Parameter( QStringLiteral( "point" ) ), fcnRasterValue, QStringLiteral( "Rasters" ) )
         << new QgsStaticExpressionFunction( QStringLiteral( "raster_attributes" ), QgsExpressionFunction::ParameterList() << QgsExpressionFunction::Parameter( QStringLiteral( "layer" ) ) << QgsExpressionFunction::Parameter( QStringLiteral( "band" ) ) << QgsExpressionFunction::Parameter( QStringLiteral( "point" ) ), fcnRasterAttributes, QStringLiteral( "Rasters" ) )
 
-        // mesh
-        << new QgsStaticExpressionFunction( QStringLiteral( "mesh_contour" ), QgsExpressionFunction::ParameterList() << QgsExpressionFunction::Parameter( QStringLiteral( "point" ), false ) << QgsExpressionFunction::Parameter( QStringLiteral( "timestamp" ), true ) << QgsExpressionFunction::Parameter( QStringLiteral( "layer" ), true ), fcnMeshContour, QStringLiteral( "Meshes" ) )
-        << new QgsStaticExpressionFunction( QStringLiteral( "mesh_data" ), QgsExpressionFunction::ParameterList() << QgsExpressionFunction::Parameter( QStringLiteral( "point" ), false ) << QgsExpressionFunction::Parameter( QStringLiteral( "dataset_group" ), false ) << QgsExpressionFunction::Parameter( QStringLiteral( "timestamp" ), true ) << QgsExpressionFunction::Parameter( QStringLiteral( "layer" ), true ), fcnMeshData, QStringLiteral( "Meshes" ) )
-
         // functions for arrays
         << new QgsArrayForeachExpressionFunction()
         << new QgsArrayFilterExpressionFunction()
diff --git a/src/core/maprenderer/qgsmaprendererjob.cpp b/src/core/maprenderer/qgsmaprendererjob.cpp
index b515aeaba92..e9057357bf0 100644
--- a/src/core/maprenderer/qgsmaprendererjob.cpp
+++ b/src/core/maprenderer/qgsmaprendererjob.cpp
@@ -47,6 +47,8 @@
 #include "qgssettingsentryimpl.h"
 #include "qgssettingstree.h"
 #include "qgsruntimeprofiler.h"
+#include "qgsmeshlayer.h"
+#include "qgsmeshlayerlabeling.h"
 
 const QgsSettingsEntryBool *QgsMapRendererJob::settingsLogCanvasRefreshEvent = new QgsSettingsEntryBool( QStringLiteral( "logCanvasRefreshEvent" ), QgsSettingsTree::sTreeMap, false );
 
@@ -253,6 +255,16 @@ bool QgsMapRendererJob::prepareLabelCache() const
         break;
       }
 
+      case Qgis::LayerType::Mesh:
+      {
+        QgsMeshLayer *l = qobject_cast< QgsMeshLayer *>( ml );
+        if ( l->labelsEnabled() && l->labeling()->requiresAdvancedEffects() )
+        {
+          canCache = false;
+        }
+        break;
+      }
+
       case Qgis::LayerType::VectorTile:
       {
         // TODO -- add detection of advanced labeling effects for vector tile layers
@@ -262,7 +274,6 @@ bool QgsMapRendererJob::prepareLabelCache() const
       case Qgis::LayerType::Raster:
       case Qgis::LayerType::Annotation:
       case Qgis::LayerType::Plugin:
-      case Qgis::LayerType::Mesh:
       case Qgis::LayerType::PointCloud:
       case Qgis::LayerType::Group:
       case Qgis::LayerType::TiledScene:
diff --git a/src/core/mesh/qgsmeshdataprovider.h b/src/core/mesh/qgsmeshdataprovider.h
index 10ecfceffeb..78d442c3908 100644
--- a/src/core/mesh/qgsmeshdataprovider.h
+++ b/src/core/mesh/qgsmeshdataprovider.h
@@ -114,6 +114,9 @@ struct CORE_EXPORT QgsMesh
   QVector<QgsMeshFace> faces SIP_SKIP;
 };
 
+Q_DECLARE_METATYPE( QgsMesh * );
+
+
 /**
  * \ingroup core
  *
diff --git a/src/core/mesh/qgsmeshlayer.cpp b/src/core/mesh/qgsmeshlayer.cpp
index b34b64a4f49..c6d732b0426 100644
--- a/src/core/mesh/qgsmeshlayer.cpp
+++ b/src/core/mesh/qgsmeshlayer.cpp
@@ -140,6 +140,12 @@ QgsMeshLayer *QgsMeshLayer::clone() const
   layer->mElevationProperties = mElevationProperties->clone();
   layer->mElevationProperties->setParent( layer );
 
+  if ( auto *lLabeling = labeling() )
+  {
+    layer->setLabeling( lLabeling->clone() );
+  }
+  layer->setLabelsEnabled( labelsEnabled() );
+
   return layer;
 }
 
@@ -1727,6 +1733,20 @@ bool QgsMeshLayer::readSymbology( const QDomNode &node, QString &errorMessage,
     setBlendMode( QgsPainting::getCompositionMode( static_cast< Qgis::BlendMode >( e.text().toInt() ) ) );
   }
 
+  // read labeling definition
+  if ( categories.testFlag( Labeling ) )
+  {
+    QgsReadWriteContextCategoryPopper p = context.enterCategory( tr( "Labeling" ) );
+
+    QDomElement labelingElement = node.firstChildElement( QStringLiteral( "labeling" ) );
+    if ( !labelingElement.isNull() )
+    {
+      QgsAbstractMeshLayerLabeling *labeling = QgsAbstractMeshLayerLabeling::create( labelingElement, context );
+      mLabelsEnabled = node.toElement().attribute( QStringLiteral( "labelsEnabled" ), QStringLiteral( "0" ) ).toInt();
+      setLabeling( labeling );
+    }
+  }
+
   // get and set the layer transparency
   if ( categories.testFlag( Rendering ) )
   {
@@ -2184,4 +2204,5 @@ void QgsMeshLayer::setLabeling( QgsAbstractMeshLayerLabeling *labeling )
 
   delete mLabeling;
   mLabeling = labeling;
+  triggerRepaint();
 }
diff --git a/src/core/mesh/qgsmeshlayer.h b/src/core/mesh/qgsmeshlayer.h
index c6661c542cc..75b002b5f1d 100644
--- a/src/core/mesh/qgsmeshlayer.h
+++ b/src/core/mesh/qgsmeshlayer.h
@@ -918,7 +918,7 @@ class CORE_EXPORT QgsMeshLayer : public QgsMapLayer, public QgsAbstractProfileSo
      * Access to labeling configuration. May be NULLPTR if labeling is not used.
      * \note Labels will only be rendered if labelsEnabled() returns TRUE.
      * \see labelsEnabled()
-     * \since QGIS 3.26
+     * \since QGIS 3.36
      */
     QgsAbstractMeshLayerLabeling *labeling() { return mLabeling; }
 
diff --git a/src/core/mesh/qgsmeshlayerlabeling.cpp b/src/core/mesh/qgsmeshlayerlabeling.cpp
index 4d178274b7a..5a8cab7d9ea 100644
--- a/src/core/mesh/qgsmeshlayerlabeling.cpp
+++ b/src/core/mesh/qgsmeshlayerlabeling.cpp
@@ -51,8 +51,9 @@ QgsPalLayerSettings QgsAbstractMeshLayerLabeling::defaultSettingsForLayer( const
 ///
 
 
-QgsMeshLayerSimpleLabeling::QgsMeshLayerSimpleLabeling( const QgsPalLayerSettings &settings )
+QgsMeshLayerSimpleLabeling::QgsMeshLayerSimpleLabeling( const QgsPalLayerSettings &settings, bool labelFaces )
   : mSettings( new QgsPalLayerSettings( settings ) )
+  , mLabelFaces( labelFaces )
 {
 }
 
@@ -63,18 +64,19 @@ QString QgsMeshLayerSimpleLabeling::type() const
 
 QgsMeshLayerSimpleLabeling *QgsMeshLayerSimpleLabeling::clone() const
 {
-  return new QgsMeshLayerSimpleLabeling( *mSettings );
+  return new QgsMeshLayerSimpleLabeling( *mSettings, mLabelFaces );
 }
 
 QgsMeshLayerLabelProvider *QgsMeshLayerSimpleLabeling::provider( QgsMeshLayer *layer ) const
 {
-  return new QgsMeshLayerLabelProvider( layer, QString(), mSettings.get() );
+  return new QgsMeshLayerLabelProvider( layer, QString(), mSettings.get(), QString(), mLabelFaces );
 }
 
 QDomElement QgsMeshLayerSimpleLabeling::save( QDomDocument &doc, const QgsReadWriteContext &context ) const
 {
   QDomElement elem = doc.createElement( QStringLiteral( "labeling" ) );
   elem.setAttribute( QStringLiteral( "type" ), QStringLiteral( "simple" ) );
+  elem.setAttribute( QStringLiteral( "labelFaces" ), mLabelFaces ? QStringLiteral( "1" ) : QStringLiteral( "0" ) );
   elem.appendChild( mSettings->writeXml( doc, context ) );
   return elem;
 }
@@ -108,7 +110,8 @@ QgsMeshLayerSimpleLabeling *QgsMeshLayerSimpleLabeling::create( const QDomElemen
   {
     QgsPalLayerSettings settings;
     settings.readXml( settingsElem, context );
-    return new QgsMeshLayerSimpleLabeling( settings );
+    const bool labelFaces = element.attribute( QStringLiteral( "labelFaces" ), QStringLiteral( "0" ) ).toInt();
+    return new QgsMeshLayerSimpleLabeling( settings, labelFaces );
   }
 
   return new QgsMeshLayerSimpleLabeling( QgsPalLayerSettings() );
diff --git a/src/core/mesh/qgsmeshlayerlabeling.h b/src/core/mesh/qgsmeshlayerlabeling.h
index c8cc2100f6f..6c649fc7a72 100644
--- a/src/core/mesh/qgsmeshlayerlabeling.h
+++ b/src/core/mesh/qgsmeshlayerlabeling.h
@@ -39,7 +39,7 @@ class QgsStyleEntityVisitorInterface;
  * \ingroup core
  * \brief Abstract base class - its implementations define different approaches to the labeling of a mesh layer.
  *
- * \since QGIS 3.38
+ * \since QGIS 3.36
  */
 class CORE_EXPORT QgsAbstractMeshLayerLabeling
 {
@@ -78,8 +78,6 @@ class CORE_EXPORT QgsAbstractMeshLayerLabeling
      *
      * \param settings Pal layer settings
      * \param providerId The id of the provider
-     *
-     * \since QGIS 3.0
      */
     virtual void setSettings( QgsPalLayerSettings *settings SIP_TRANSFER, const QString &providerId = QString() ) = 0;
 
@@ -95,8 +93,6 @@ class CORE_EXPORT QgsAbstractMeshLayerLabeling
      *
      * This method multiplies the opacity of the labeling elements (text, shadow, buffer etc.)
      * by \a opacity effectively changing the opacity of the whole labeling elements.
-     *
-     * \since QGIS 3.32
      */
     virtual void multiplyOpacity( double opacityFactor ) { Q_UNUSED( opacityFactor ); };
 
@@ -144,13 +140,17 @@ class CORE_EXPORT QgsAbstractMeshLayerLabeling
  * \ingroup core
  * \brief Basic implementation of the labeling interface for mesh layer.
  *
- * \since QGIS 3.38
+ * \since QGIS 3.36
  */
 class CORE_EXPORT QgsMeshLayerSimpleLabeling : public QgsAbstractMeshLayerLabeling
 {
   public:
-    //! Constructs simple labeling configuration with given initial settings
-    explicit QgsMeshLayerSimpleLabeling( const QgsPalLayerSettings &settings );
+
+    /**
+     *  Constructs simple labeling configuration with given initial \a settings.
+     *  Labels are placed on mesh vertices unless \a labelFaces is TRUE, when they are placed on mesh faces.
+     */
+    explicit QgsMeshLayerSimpleLabeling( const QgsPalLayerSettings &settings, bool labelFaces = false );
 
     QString type() const override;
     QgsMeshLayerSimpleLabeling *clone() const override SIP_FACTORY;
@@ -175,6 +175,7 @@ class CORE_EXPORT QgsMeshLayerSimpleLabeling : public QgsAbstractMeshLayerLabeli
 
   private:
     std::unique_ptr<QgsPalLayerSettings> mSettings;
+    bool mLabelFaces = false;
 };
 
 #endif // QGSMESHLAYERLABELING_H
diff --git a/src/core/mesh/qgsmeshlayerlabelprovider.cpp b/src/core/mesh/qgsmeshlayerlabelprovider.cpp
index fa09180e277..c08289d5bcb 100644
--- a/src/core/mesh/qgsmeshlayerlabelprovider.cpp
+++ b/src/core/mesh/qgsmeshlayerlabelprovider.cpp
@@ -49,9 +49,10 @@
 
 using namespace pal;
 
-QgsMeshLayerLabelProvider::QgsMeshLayerLabelProvider( QgsMeshLayer *layer, const QString &providerId, const QgsPalLayerSettings *settings, const QString &layerName )
+QgsMeshLayerLabelProvider::QgsMeshLayerLabelProvider( QgsMeshLayer *layer, const QString &providerId, const QgsPalLayerSettings *settings, const QString &layerName, bool labelFaces )
   : QgsAbstractLabelProvider( layer, providerId )
   , mSettings( settings ? * settings : QgsPalLayerSettings() )
+  , mLabelFaces( labelFaces )
   , mCrs( layer->crs() )
 {
   mName = layerName.isEmpty() ? layer->id() : layerName;
@@ -73,38 +74,13 @@ void QgsMeshLayerLabelProvider::init()
 
   mPriority = 1 - mSettings.priority / 10.0; // convert 0..10 --> 1..0
 
-  if ( mLabelFaces )
-  {
-    mVectorLayer = std::make_unique<QgsVectorLayer>( QStringLiteral( "Polygon?crs=%1" ).arg( mCrs.authid() ), QStringLiteral( "faces" ), QStringLiteral( "memory" ) );
-    QgsMesh *mesh = qobject_cast<QgsMeshLayer *>( layer() )->nativeMesh();
-    int faceCount = mesh->faceCount();
-    QList<QgsFeature> features;
-    for ( int i = 0; i < faceCount; i++ )
-    {
-      QgsFeature f;
-      QgsGeometry geom = QgsMeshUtils::toGeometry( mesh->face( i ), mesh->vertices );
-      f.setGeometry( geom );
-      features << f;
-    }
-    mVectorLayer->dataProvider()->addFeatures( features );
-  }
-  else
-  {
-    mVectorLayer = std::make_unique<QgsVectorLayer>( QStringLiteral( "Polygon?crs=%1" ).arg( mCrs.authid() ), QStringLiteral( "faces" ), QStringLiteral( "memory" ) );
-    QgsMesh *mesh = qobject_cast<QgsMeshLayer *>( layer() )->nativeMesh();
-    int vertexCount = mesh->vertexCount();
-    QList<QgsFeature> features;
-    for ( int i = 0; i < vertexCount; i++ )
-    {
-      QgsFeature f;
-      QgsGeometry geom = QgsGeometry( new QgsPoint( mesh->vertex( i ) ) );
-      f.setGeometry( geom );
-      features << f;
-    }
-    mVectorLayer->dataProvider()->addFeatures( features );
-  }
-
-  mVectorLabelProvider = std::make_unique<QgsVectorLayerLabelProvider>( mVectorLayer.get(), QStringLiteral(), false, &mSettings );
+  mVectorLabelProvider = std::make_unique<QgsVectorLayerLabelProvider>(
+                           mLabelFaces ? Qgis::GeometryType::Polygon : Qgis::GeometryType::Point,
+                           QgsFields(),
+                           mCrs,
+                           QString(),
+                           &mSettings,
+                           mLayer );
 
   if ( mLabelFaces )
   {
@@ -129,6 +105,8 @@ bool QgsMeshLayerLabelProvider::prepare( QgsRenderContext &context, QSet<QString
 {
   const QgsMapSettings &mapSettings = mEngine->mapSettings();
 
+  mVectorLabelProvider->setEngine( mEngine );
+
   return mSettings.prepare( context, attributeNames, QgsFields(), mapSettings, mCrs );
 }
 
diff --git a/src/core/mesh/qgsmeshlayerlabelprovider.h b/src/core/mesh/qgsmeshlayerlabelprovider.h
index 3bd4d13053e..b91a0f8a3fb 100644
--- a/src/core/mesh/qgsmeshlayerlabelprovider.h
+++ b/src/core/mesh/qgsmeshlayerlabelprovider.h
@@ -35,7 +35,7 @@ class QgsMeshLayer;
  *
  * \note this class is not a part of public API yet. See notes in QgsLabelingEngine
  * \note not available in Python bindings
- * \since QGIS 2.12
+ * \since QGIS 3.36
  */
 class CORE_EXPORT QgsMeshLayerLabelProvider : public QgsAbstractLabelProvider
 {
@@ -45,7 +45,8 @@ class CORE_EXPORT QgsMeshLayerLabelProvider : public QgsAbstractLabelProvider
     explicit QgsMeshLayerLabelProvider( QgsMeshLayer *layer,
                                         const QString &providerId,
                                         const QgsPalLayerSettings *settings,
-                                        const QString &layerName = QString() );
+                                        const QString &layerName = QString(),
+                                        bool labelFaces = false );
 
     ~QgsMeshLayerLabelProvider() override;
 
@@ -86,6 +87,9 @@ class CORE_EXPORT QgsMeshLayerLabelProvider : public QgsAbstractLabelProvider
      */
     const QgsPalLayerSettings &settings() const;
 
+    //! Returns FALSE if labeling mesh vertices, TRUE if labeling mesh faces
+    bool labelFaces() const { return mLabelFaces; }
+
   protected:
     //! initialization method - called from constructors
     void init();
@@ -106,7 +110,6 @@ class CORE_EXPORT QgsMeshLayerLabelProvider : public QgsAbstractLabelProvider
 
   private:
     std::unique_ptr<QgsVectorLayerLabelProvider> mVectorLabelProvider;
-    std::unique_ptr<QgsVectorLayer> mVectorLayer;
 
     friend class TestQgsLabelingEngine;
 };
diff --git a/src/core/mesh/qgsmeshlayerrenderer.cpp b/src/core/mesh/qgsmeshlayerrenderer.cpp
index 3e84998364d..8ea620ff673 100644
--- a/src/core/mesh/qgsmeshlayerrenderer.cpp
+++ b/src/core/mesh/qgsmeshlayerrenderer.cpp
@@ -34,11 +34,14 @@
 #include "qgsmeshlayerinterpolator.h"
 #include "qgsmeshlayerutils.h"
 #include "qgsmeshvectorrenderer.h"
+#include "qgsmeshlayerlabeling.h"
+#include "qgsmeshlayerlabelprovider.h"
 #include "qgsmapclippingutils.h"
 #include "qgscolorrampshader.h"
 #include "qgsmaplayerelevationproperties.h"
 #include "qgsapplication.h"
 #include "qgsruntimeprofiler.h"
+#include "qgsexpressioncontextutils.h"
 
 QgsMeshLayerRenderer::QgsMeshLayerRenderer(
   QgsMeshLayer *layer,
@@ -79,6 +82,9 @@ QgsMeshLayerRenderer::QgsMeshLayerRenderer(
 
   calculateOutputSize();
 
+  QSet<QString> attrs;
+  prepareLabeling( layer, attrs );
+
   mClippingRegions = QgsMapClippingUtils::collectClippingRegionsForLayer( *renderContext(), layer );
 
   if ( layer->elevationProperties() && layer->elevationProperties()->hasElevation() )
@@ -333,6 +339,8 @@ bool QgsMeshLayerRenderer::render()
   renderMesh();
   renderVectorDataset();
 
+  registerLabelFeatures();
+
   return !renderContext()->renderingStopped();
 }
 
@@ -648,3 +656,89 @@ void QgsMeshLayerRenderer::renderVectorDataset()
     renderer->draw();
 }
 
+void QgsMeshLayerRenderer::prepareLabeling( QgsMeshLayer *layer, QSet<QString> &attributeNames )
+{
+  QgsRenderContext &context = *renderContext();
+
+  if ( QgsLabelingEngine *engine = context.labelingEngine() )
+  {
+    if ( layer->labelsEnabled() )
+    {
+      mLabelProvider = layer->labeling()->provider( layer );
+      if ( mLabelProvider )
+      {
+        auto c = context.expressionContext();
+
+        c.appendScope( QgsExpressionContextUtils::meshExpressionScope( mLabelProvider->labelFaces() ? QgsMesh::Face : QgsMesh::Vertex ) );
+        c.lastScope()->setVariable( QStringLiteral( "_native_mesh" ), QVariant::fromValue( &mNativeMesh ) );
+        context.setExpressionContext( c );
+
+        engine->addProvider( mLabelProvider );
+        if ( !mLabelProvider->prepare( context, attributeNames ) )
+        {
+          engine->removeProvider( mLabelProvider );
+          mLabelProvider = nullptr; // deleted by engine
+        }
+      }
+    }
+  }
+}
+
+void QgsMeshLayerRenderer::registerLabelFeatures()
+{
+  if ( !mLabelProvider )
+    return;
+
+  QgsRenderContext &context = *renderContext();
+
+  QgsExpressionContextScope *scope = context.expressionContext().activeScopeForVariable( QStringLiteral( "_native_mesh" ) );
+
+  const QList<int> trianglesInExtent = mTriangularMesh.faceIndexesForRectangle( renderContext()->mapExtent() );
+
+  if ( mLabelProvider->labelFaces() )
+  {
+    if ( !mTriangularMesh.contains( QgsMesh::ElementType::Face ) )
+      return;
+    const QSet<int> nativeFacesInExtent = QgsMeshUtils::nativeFacesFromTriangles( trianglesInExtent,
+                                          mTriangularMesh.trianglesToNativeFaces() );
+
+    for ( const int i : nativeFacesInExtent )
+    {
+      if ( context.renderingStopped() )
+        break;
+
+      if ( i < 0 || i >= mNativeMesh.faces.count() )
+        continue;
+
+      scope->setVariable( QStringLiteral( "_mesh_face_index" ), i, false );
+
+      QgsFeature f( i );
+      QgsGeometry geom = QgsMeshUtils::toGeometry( mNativeMesh.face( i ), mNativeMesh.vertices );
+      f.setGeometry( geom );
+      mLabelProvider->registerFeature( f, context );
+    }
+  }
+  else
+  {
+    if ( !mTriangularMesh.contains( QgsMesh::ElementType::Vertex ) )
+      return;
+    const QSet<int> nativeVerticesInExtent = QgsMeshUtils::nativeVerticesFromTriangles( trianglesInExtent,
+        mTriangularMesh.triangles() );
+
+    for ( const int i : nativeVerticesInExtent )
+    {
+      if ( context.renderingStopped() )
+        break;
+
+      if ( i < 0 || i >= mNativeMesh.vertexCount() )
+        continue;
+
+      scope->setVariable( QStringLiteral( "_mesh_vertex_index" ), i, false );
+
+      QgsFeature f( i );
+      QgsGeometry geom = QgsGeometry( new QgsPoint( mNativeMesh.vertex( i ) ) );
+      f.setGeometry( geom );
+      mLabelProvider->registerFeature( f, context );
+    }
+  }
+}
diff --git a/src/core/mesh/qgsmeshlayerrenderer.h b/src/core/mesh/qgsmeshlayerrenderer.h
index aa9023e9ade..e7e9100d512 100644
--- a/src/core/mesh/qgsmeshlayerrenderer.h
+++ b/src/core/mesh/qgsmeshlayerrenderer.h
@@ -36,6 +36,7 @@ class QgsMeshLayer;
 #include "qgsmapclippingregion.h"
 
 class QgsRenderContext;
+class QgsMeshLayerLabelProvider;
 
 ///@cond PRIVATE
 
@@ -100,12 +101,14 @@ class QgsMeshLayerRenderer : public QgsMapLayerRenderer
     bool forceRasterRender() const override;
 
   private:
+    void prepareLabeling( QgsMeshLayer *layer, QSet<QString> &attributeNames );
     void renderMesh();
     void renderEdgeMesh( const QgsMeshRendererMeshSettings &settings, const QList<int> &edgesInExtent );
     void renderFaceMesh( const QgsMeshRendererMeshSettings &settings, const QVector<QgsMeshFace> &faces, const QList<int> &facesInExtent );
     void renderScalarDataset();
     void renderScalarDatasetOnEdges( const QgsMeshRendererScalarSettings &scalarSettings );
     void renderScalarDatasetOnFaces( const QgsMeshRendererScalarSettings &scalarSettings );
+    void registerLabelFeatures();
 
     void renderVectorDataset();
     void copyTriangularMeshes( QgsMeshLayer *layer, QgsRenderContext &context );
@@ -117,6 +120,12 @@ class QgsMeshLayerRenderer : public QgsMapLayerRenderer
     QColor colorAt( QgsColorRampShader *shader, double val ) const;
     bool mIsEditable = false;
 
+    /**
+     * used with new labeling engine (QgsLabelingEngine): provider for labels.
+     * may be NULLPTR. no need to delete: if exists it is owned by labeling engine
+     */
+    QgsMeshLayerLabelProvider *mLabelProvider = nullptr;
+
   protected:
     QString mLayerName;
 
diff --git a/tests/src/core/testqgsexpression.cpp b/tests/src/core/testqgsexpression.cpp
index 81f20fe9905..54185e6f1b1 100644
--- a/tests/src/core/testqgsexpression.cpp
+++ b/tests/src/core/testqgsexpression.cpp
@@ -2016,13 +2016,6 @@ class TestQgsExpression: public QObject
       QTest::newRow( "raster_attributes wrong band 2" ) << QStringLiteral( "raster_attributes('%1',2,243)" ).arg( mRasterLayerWithAttributeTable->name() ) << true << QVariant();
       QTest::newRow( "raster_attributes no attributes" ) << QStringLiteral( "raster_attributes('%1',1,1)" ).arg( mRasterLayerWithAttributeTable->name() ) << false << QVariant();
 
-      // mesh_contour tests
-      QTest::newRow( "mesh_contour invalid geometry" ) << QStringLiteral( "mesh_contour('invalid geom', make_datetime(2023,11,29,10,34,52), '%1')" ).arg( mMeshLayer->name() ) << true << QVariant();
-      QTest::newRow( "mesh_contour valid" ) << QStringLiteral( "mesh_contour(make_point(2000,3000), make_datetime(2023,11,29,10,34,52), '%1')" ).arg( mMeshLayer->name() ) << false << QVariant( 200.0 );
-      QTest::newRow( "mesh_data invalid geometry" ) << QStringLiteral( "mesh_data('invalid geom', 'Bed Elevation', make_datetime(2023,11,29,10,34,52), '%1')" ).arg( mMeshLayer->name() ) << true << QVariant();
-      QTest::newRow( "mesh_data invalid group" ) << QStringLiteral( "mesh_data(make_point(2000,3000), 'Bad elevation', make_datetime(2023,11,29,10,34,52), '%1')" ).arg( mMeshLayer->name() ) << true << QVariant();
-      QTest::newRow( "mesh_data valid" ) << QStringLiteral( "mesh_data(make_point(2000,3000), 'Bed Elevation', make_datetime(2023,11,29,10,34,52), '%1')" ).arg( mMeshLayer->name() ) << false << QVariant( 200.0 );
-
       //test conversions to bool
       QTest::newRow( "feature to bool false" ) << QStringLiteral( "case when get_feature('none','none',499) then true else false end" ) << false << QVariant( false );
       QTest::newRow( "feature to bool true" ) << QStringLiteral( "case when get_feature('test','col1',10) then true else false end" ) << false << QVariant( true );