diff --git a/external/mdal/frmts/mdal_mike21.cpp b/external/mdal/frmts/mdal_mike21.cpp
new file mode 100644
index 00000000000..1f63758d62e
--- /dev/null
+++ b/external/mdal/frmts/mdal_mike21.cpp
@@ -0,0 +1,491 @@
+/*
+ MDAL - Mesh Data Abstraction Library (MIT License)
+ Copyright (C) 2023 Lutra Consulting Ltd.
+*/
+
+#include <stddef.h>
+#include <iosfwd>
+#include <iostream>
+#include <fstream>
+#include <sstream>
+#include <string>
+#include <vector>
+#include <map>
+#include <cassert>
+#include <limits>
+#include <algorithm>
+#include <regex>
+
+#include "mdal_mike21.hpp"
+#include "mdal.h"
+#include "mdal_utils.hpp"
+#include "mdal_logger.hpp"
+
+#define DRIVER_NAME "Mike21"
+
+// function to split using regex, by default split on whitespace characters
+std::vector<std::string> regex_split( const std::string &input, const std::regex &split_regex = std::regex{"\\s+"} )
+{
+  std::sregex_token_iterator iter( input.begin(), input.end(), split_regex, -1 );
+  std::sregex_token_iterator end;
+  return {iter, end};
+}
+
+static bool parse_vertex_id_gaps( std::map<size_t, size_t> &vertexIDtoIndex, size_t vertexIndex, size_t vertexID )
+{
+  if ( vertexIndex == vertexID )
+    return false;
+
+  std::map<size_t, size_t>::iterator search = vertexIDtoIndex.find( vertexID );
+  if ( search != vertexIDtoIndex.end() )
+  {
+    MDAL::Log::warning( Warn_ElementNotUnique, DRIVER_NAME, "could not find vertex" );
+    return true;
+  }
+
+  vertexIDtoIndex[vertexID] = vertexIndex;
+  return false;
+}
+
+static void persist_native_index( std::vector<double> &arr, size_t nativeID, size_t ourId, size_t maxOurId )
+{
+  if ( !arr.empty() || ( nativeID != ourId + 1 ) )
+  {
+    // we have gaps in face indexing
+    if ( arr.empty() )
+    {
+      arr.resize( maxOurId );
+      for ( size_t i = 0; i < ourId; ++i )
+        arr[i] = static_cast<double>( i + 1 );
+    }
+    arr[ourId] = static_cast<double>( nativeID );
+  }
+}
+
+MDAL::MeshMike21::MeshMike21( size_t faceVerticesMaximumCount,
+                              const std::string &uri,
+                              const std::map<size_t, size_t> vertexIDtoIndex )
+  : MemoryMesh( DRIVER_NAME,
+                faceVerticesMaximumCount,
+                uri )
+  , mVertexIDtoIndex( vertexIDtoIndex )
+{
+}
+
+MDAL::MeshMike21::~MeshMike21() = default;
+
+size_t MDAL::MeshMike21::vertexIndex( size_t vertexID ) const
+{
+  auto ni2i = mVertexIDtoIndex.find( vertexID );
+  if ( ni2i != mVertexIDtoIndex.end() )
+  {
+    return  ni2i->second; // convert from ID to index
+  }
+  return vertexID;
+}
+
+size_t MDAL::MeshMike21::maximumVertexId() const
+{
+  size_t maxIndex = verticesCount() - 1;
+  if ( mVertexIDtoIndex.empty() )
+    return maxIndex;
+  else
+  {
+    // std::map is sorted!
+    size_t maxID = mVertexIDtoIndex.rbegin()->first;
+    return std::max( maxIndex, maxID );
+  }
+}
+
+MDAL::DriverMike21::DriverMike21( ):
+  Driver( DRIVER_NAME,
+          "Mike21 Mesh File",
+          "*.mesh",
+          Capability::ReadMesh | Capability::SaveMesh
+        )
+{
+}
+
+MDAL::DriverMike21 *MDAL::DriverMike21::create()
+{
+  return new DriverMike21();
+}
+
+MDAL::DriverMike21::~DriverMike21() = default;
+
+bool MDAL::DriverMike21::canReadHeader( const std::string &line )
+{
+  bool header2012 = std::regex_match( line, mRegexHeader2012 );
+  bool header2011 = std::regex_match( line, mRegexHeader2011 );
+  return header2011 || header2012;
+}
+
+bool MDAL::DriverMike21::canReadMesh( const std::string &uri )
+{
+  std::ifstream in = MDAL::openInputFile( uri );
+  std::string line;
+
+  if ( !MDAL::getHeaderLine( in, line ) || !canReadHeader( line ) ||  !MDAL::contains( filters(), MDAL::fileExtension( uri ) ) )
+  {
+    return false;
+  }
+
+  return true;
+}
+
+size_t MDAL::DriverMike21::getVertexCount( const std::string &line )
+{
+  auto matchResults = std::smatch{};
+  if ( std::regex_search( line, matchResults, mRegexHeader2012 ) )
+  {
+    if ( matchResults.size() > 4 )
+    {
+      return std::stoi( matchResults[3].str() );
+    }
+  }
+
+  if ( std::regex_search( line, matchResults, mRegexHeader2011 ) )
+  {
+    if ( matchResults.size() > 2 )
+    {
+      return std::stoi( matchResults[1].str() );
+    }
+  }
+
+  return 0;
+}
+
+std::string MDAL::DriverMike21::getCrs( const std::string &line )
+{
+  auto matchResults = std::smatch{};
+  if ( std::regex_search( line, matchResults, mRegexHeader2012 ) )
+  {
+    if ( matchResults.size() > 5 )
+    {
+      return matchResults[4].str();
+    }
+  }
+
+  if ( std::regex_search( line, matchResults, mRegexHeader2011 ) )
+  {
+    if ( matchResults.size() > 3 )
+    {
+      return matchResults[2].str();
+    }
+  }
+
+  return "";
+}
+
+std::unique_ptr<MDAL::Mesh> MDAL::DriverMike21::load( const std::string &meshFile, const std::string & )
+{
+  mMeshFile = meshFile;
+
+  MDAL::Log::resetLastStatus();
+
+  std::ifstream in = MDAL::openInputFile( meshFile );
+
+  std::string line;
+  if ( !std::getline( in, line ) || !canReadHeader( line ) )
+  {
+    MDAL::Log::error( MDAL_Status::Err_UnknownFormat, name(), meshFile + " could not be opened" );
+    return nullptr;
+  }
+
+  std::string crs = getCrs( line );
+
+  size_t vertexCount = getVertexCount( line );
+  size_t faceCount = 0;
+  size_t maxVerticesPerFace = 2;
+
+  size_t lineNumber = 1;
+
+  while ( std::getline( in, line ) )
+  {
+    if ( lineNumber == vertexCount + 1 )
+    {
+      auto matchResults = std::smatch{};
+      if ( std::regex_search( line, matchResults, mRegexElementHeader ) )
+      {
+        if ( matchResults.size() >= 4 )
+        {
+          faceCount = MDAL::toSizeT( matchResults[1].str() );
+          maxVerticesPerFace = MDAL::toSizeT( matchResults[2].str() );
+          size_t meshType = MDAL::toSizeT( matchResults[3].str() );
+
+          if ( !( meshType == 21 || meshType == 25 ) )
+          {
+            MDAL::Log::error( MDAL_Status::Err_InvalidData, name(), "unknow mesh type." );
+            return nullptr;
+          }
+        }
+        else
+        {
+          MDAL::Log::error( MDAL_Status::Err_InvalidData, name(), "element header not in valid format." );
+          return nullptr;
+        }
+      }
+      else
+      {
+        MDAL::Log::error( MDAL_Status::Err_InvalidData, name(), "element header not in valid format." );
+        return nullptr;
+      }
+
+    }
+    lineNumber++;
+  }
+
+  // number of lines in file does not match number of vertices and faces specifed in first and element line
+  if ( lineNumber > 2 + vertexCount + faceCount )
+  {
+    MDAL::Log::error( MDAL_Status::Err_InvalidData, name(), "Number of lines in file does not fit with number of vertexes and faces specified." );
+    return nullptr;
+  }
+
+  in.clear();
+  in.seekg( 0, std::ios::beg );
+
+  Vertices vertices( vertexCount );
+  Faces faces( faceCount );
+
+  std::map<size_t, size_t> vertexIDtoIndex;
+  std::vector<double> vertexType( vertexCount );
+
+  std::vector<double> nativeVertexIds;
+  std::vector<double> nativeFaceIds;
+
+  size_t lastVertexID = 0;
+  size_t faceIndex = 0;
+  size_t vertexIndex = 0;
+
+  std::vector<std::string> chunks;
+  lineNumber = 0;
+
+  while ( std::getline( in, line ) )
+  {
+    if ( 0 < lineNumber && lineNumber < vertexCount + 1 )
+    {
+      chunks = regex_split( MDAL::trim( line ) );
+      if ( chunks.size() != 5 )
+      {
+        MDAL::Log::error( MDAL_Status::Err_InvalidData, name(), "vertex line in invalid format." );
+        return nullptr;
+      }
+
+      size_t nodeID = toSizeT( chunks[0] );
+
+      if ( nodeID != 0 )
+      {
+        // specification of Mike21 does not state if vertexIDs need to continuos, expect that they might be not
+        // in the same way as in 2DM
+        if ( ( lastVertexID != 0 ) && ( nodeID <= lastVertexID ) )
+        {
+          // the algorithm requires that the file has points orderer by index
+          MDAL::Log::error( MDAL_Status::Err_InvalidData, name(), "nodes are not ordered by index" );
+          return nullptr;
+        }
+        lastVertexID = nodeID;
+      }
+
+      // in case we have gaps/reorders in native indexes, store it
+      persist_native_index( nativeVertexIds, nodeID, vertexIndex, vertexCount );
+      parse_vertex_id_gaps( vertexIDtoIndex, vertexIndex, nodeID - 1 );
+
+      assert( vertexIndex < vertexCount );
+      Vertex &vertex = vertices[vertexIndex];
+      vertex.x = toDouble( chunks[1] );
+      vertex.y = toDouble( chunks[2] );
+      vertex.z = toDouble( chunks[3] );
+      vertexType[vertexIndex] = MDAL::toInt( chunks[4] );
+      vertexIndex++;
+    }
+
+    if ( vertexCount + 1 < lineNumber )
+    {
+      chunks = regex_split( MDAL::trim( line ) );
+      assert( faceIndex < faceCount );
+
+      size_t faceVertexCount = chunks.size() - 1;
+      // if the face should have 4 vertexes last chunk has value 0
+      // it actually means that there are only 3 vertexes
+      if ( faceVertexCount == 4 && chunks.size() == 5 )
+      {
+        if ( MDAL::toSizeT( chunks[4] ) == 0 )
+          faceVertexCount = faceVertexCount - 1;
+      }
+
+      assert( ( faceVertexCount == 3 ) || ( faceVertexCount == 4 ) );
+      if ( maxVerticesPerFace < faceVertexCount )
+        maxVerticesPerFace = faceVertexCount;
+
+      Face &face = faces[faceIndex];
+      face.resize( faceVertexCount );
+
+      // in case we have gaps/reorders in native indexes, store it
+      size_t nativeID = MDAL::toSizeT( chunks[0] );
+      persist_native_index( nativeFaceIds, nativeID, faceIndex, faceCount );
+
+      for ( size_t i = 0; i < faceVertexCount; ++i )
+        face[i] = MDAL::toSizeT( chunks[i + 1] ) - 1; // Mike21 is numbered from 1
+
+      faceIndex++;
+    }
+
+    lineNumber++;
+  }
+
+  for ( std::vector<Face>::iterator it = faces.begin(); it != faces.end(); ++it )
+  {
+    Face &face = *it;
+    for ( Face::size_type nd = 0; nd < face.size(); ++nd )
+    {
+      size_t nodeID = face[nd];
+
+      std::map<size_t, size_t>::iterator ni2i = vertexIDtoIndex.find( nodeID );
+      if ( ni2i != vertexIDtoIndex.end() )
+      {
+        face[nd] = ni2i->second; // convert from ID to index
+      }
+      else if ( vertices.size() < nodeID )
+      {
+        MDAL::Log::warning( MDAL_Status::Warn_ElementWithInvalidNode, name(), "found invalid node" );
+      }
+    }
+  }
+
+  // create the mesh and set the required data
+  std::unique_ptr< MeshMike21 > mesh(
+    new MeshMike21(
+      maxVerticesPerFace,
+      mMeshFile,
+      vertexIDtoIndex
+    )
+  );
+  mesh->setFaces( std::move( faces ) );
+  mesh->setVertices( std::move( vertices ) );
+
+  // Add Vertex Type
+  MDAL::addVertexScalarDatasetGroup( mesh.get(), vertexType, "VertexType" );
+
+  // Add Bed Elevation
+  MDAL::addBedElevationDatasetGroup( mesh.get(), mesh->vertices() );
+
+  if ( !nativeFaceIds.empty() )
+    MDAL::addFaceScalarDatasetGroup( mesh.get(), nativeFaceIds, "NativeFaceIds" );
+  if ( !nativeVertexIds.empty() )
+    MDAL::addVertexScalarDatasetGroup( mesh.get(), nativeVertexIds, "NativeVertexIds" );
+
+  mesh->setSourceCrs( crs );
+
+  return std::unique_ptr<Mesh>( mesh.release() );
+}
+
+void MDAL::DriverMike21::save( const std::string &fileName, const std::string &, MDAL::Mesh *mesh )
+{
+  MDAL::Log::resetLastStatus();
+
+  std::ofstream file = MDAL::openOutputFile( fileName, std::ofstream::out );
+
+  if ( !file.is_open() )
+  {
+    MDAL::Log::error( MDAL_Status::Err_FailToWriteToDisk, name(), "Could not open file " + fileName );
+  }
+
+  std::string line = std::to_string( mesh->verticesCount() ) + " " + mesh->crs();
+  file << line << std::endl;
+
+  std::vector<double> vertexTypes;
+
+  std::shared_ptr<MDAL::DatasetGroup> vertexTypeDG = mesh->group( "VertexType" );
+  if ( vertexTypeDG )
+  {
+    vertexTypes.resize( mesh->verticesCount() );
+    auto d = vertexTypeDG->datasets[0];
+    d->scalarData( 0, mesh->verticesCount(), vertexTypes.data() );
+  }
+
+  // write vertices
+  std::unique_ptr<MDAL::MeshVertexIterator> vertexIterator = mesh->readVertices();
+  double vertex[3];
+  for ( size_t i = 0; i < mesh->verticesCount(); ++i )
+  {
+    vertexIterator->next( 1, vertex );
+    line = "";
+    line.append( std::to_string( i + 1 ) );
+    for ( size_t j = 0; j < 2; ++j )
+    {
+      line.append( " " );
+      line.append( MDAL::coordinateToString( vertex[j] ) );
+    }
+    line.append( " " );
+    line.append( MDAL::doubleToString( vertex[2] ) );
+
+    line.append( " " );
+    if ( vertexTypes.size() == mesh->verticesCount() )
+    {
+      line.append( MDAL::doubleToString( vertexTypes.at( i ) ) );
+    }
+    else
+    {
+      line.append( MDAL::doubleToString( 0 ) );
+    }
+
+    file << line << std::endl;
+  }
+
+  //write element header line
+  size_t elementType = 0;
+  if ( mesh->faceVerticesMaximumCount() == 3 )
+  {
+    elementType = 21;
+  }
+  else if ( mesh->faceVerticesMaximumCount() == 4 )
+  {
+    elementType = 25;
+  }
+
+  line = std::to_string( mesh->facesCount() );
+  line.append( " " );
+  line.append( std::to_string( mesh->faceVerticesMaximumCount() ) );
+  line.append( " " );
+  line.append( std::to_string( elementType ) );
+  file << line << std::endl;
+
+  // write faces
+  std::vector<int> vertexIndices( mesh->faceVerticesMaximumCount() );
+  std::unique_ptr<MDAL::MeshFaceIterator> faceIterator = mesh->readFaces();
+  for ( size_t i = 0; i < mesh->facesCount(); ++i )
+  {
+    int faceOffsets[1];
+    faceIterator->next( 1, faceOffsets, 4, vertexIndices.data() );
+
+    if ( faceOffsets[0] > 2 && faceOffsets[0] < 5 )
+    {
+      line = "";
+      line.append( std::to_string( i + 1 ) );
+
+      for ( int j = 0; j < faceOffsets[0]; ++j )
+      {
+        line.append( " " );
+        line.append( std::to_string( vertexIndices[j] + 1 ) );
+      }
+
+      // if face has 3 vertexes but the mesh as whole is marked as having
+      // 4 vertex at maximum, the last element should 0 - indicating no vertex there
+      if ( faceOffsets[0] == 3 && mesh->faceVerticesMaximumCount() == 4 )
+      {
+        line.append( " " );
+        line.append( "0" );
+      }
+
+    }
+    file << line << std::endl;
+  }
+
+  file.close();
+}
+
+std::string MDAL::DriverMike21::saveMeshOnFileSuffix() const
+{
+  return "mesh";
+}
diff --git a/external/mdal/frmts/mdal_mike21.hpp b/external/mdal/frmts/mdal_mike21.hpp
new file mode 100644
index 00000000000..60c36256926
--- /dev/null
+++ b/external/mdal/frmts/mdal_mike21.hpp
@@ -0,0 +1,104 @@
+/*
+ MDAL - Mesh Data Abstraction Library (MIT License)
+ Copyright (C) 2023 Lutra Consulting Ltd.
+*/
+
+#ifndef MDAL_MIKE21_HPP
+#define MDAL_MIKE21_HPP
+
+#include <string>
+#include <memory>
+#include <regex>
+
+#include "mdal_data_model.hpp"
+#include "mdal_memory_data_model.hpp"
+#include "mdal.h"
+#include "mdal_driver.hpp"
+
+namespace MDAL
+{
+  class MeshMike21: public MemoryMesh
+  {
+    public:
+      MeshMike21( size_t faceVerticesMaximumCount,
+                  const std::string &uri,
+                  const std::map<size_t, size_t> vertexIDtoIndex
+                );
+      ~MeshMike21() override;
+
+      /**
+       * Mike21 may supports gaps in the vertex indexing,
+       * but we use continuos array of vertices in MDAL
+       * \param vertexID internal index/ID of the vertex that native format uses
+       * \returns index of the vertex in the continuous array of vertices we returned by readVertices().
+       *          For invalid vertexID it is returned index that is out of vertices array bounds.
+       */
+      virtual size_t vertexIndex( size_t vertexID ) const;
+
+      /**
+       * Returns maximum vertex ID.
+       * For meshes without gaps in vertex indexing, it is vertex count - 1
+       */
+      virtual size_t maximumVertexId() const;
+
+    private:
+      /**
+       * Mike21 might supports "gaps" in the mesh indexing
+       * Store only the indices that have different index and ID
+       * https://github.com/lutraconsulting/MDAL/issues/51
+       */
+      std::map<size_t, size_t> mVertexIDtoIndex;
+  };
+
+  /**
+   * Mike21 format specification
+   * Text file format representing mesh vertices and faces
+   * Vertices are stored as - id x y z type_of_vertex (0 for water, 1 for land and above 1 for all other boundaries)
+   * The format supports triangles and quads
+   * Faces are stored as id 1 2 3 [4]
+   *
+   * full specification here: https://www.xmswiki.com/wiki/SMS:MIKE_21_*.mesh
+   *
+   * The format has two special lines in the file:
+   * First line - that can have format either:
+   *    int int int string - these are type of bathymetry data, units of bathymetry data, number of vertices, CRS
+   *    int string - these are number of vertices, CRS
+   * Line after the vertices (Element header line):
+   *    int int int - number of faces, max number of vertexes per face, mesh type code
+   *    [21 for faces with triangles only, 25 for faces with also quadrangular elements]
+   *
+   * The driver creates at least two Vertex Scalar datasets on the Mesh. "Bed Elevation" for Z coordinate
+   * and "VertexType" storing the type of vertex (mentioned above). Besides "NativeVertexIds" and "NativeFaceIds"
+   * maybe created if the IDS in input data are not continuous.
+   */
+  class DriverMike21: public Driver
+  {
+    public:
+      DriverMike21();
+      ~DriverMike21() override;
+      DriverMike21 *create() override;
+
+      int faceVerticesMaximumCount() const override {return 4;}
+
+      bool canReadMesh( const std::string &uri ) override;
+      std::unique_ptr< Mesh > load( const std::string &meshFile, const std::string &meshName = "" ) override;
+      void save( const std::string &fileName, const std::string &, Mesh *mesh ) override;
+
+      std::string saveMeshOnFileSuffix() const override;
+
+    private:
+      std::string mMeshFile;
+      // regex for header line in form of - integer string
+      const std::regex mRegexHeader2011 = std::regex( "(\\d+)\\s+(.+)(\\s+)?" );
+      // regex for header line in form of - integer integer integer string
+      const std::regex mRegexHeader2012 = std::regex( "(\\d+)\\s+(\\d+)\\s+(\\d+)\\s+(.+)(\\s+)?" );
+      // regex for element header line - integer integer integer(code with two number)
+      const std::regex mRegexElementHeader = std::regex( "(\\d+)\\s+(\\d)\\s+(\\d{2})(\\s+)?" );
+
+      bool canReadHeader( const std::string &line );
+      size_t getVertexCount( const std::string &line );
+      std::string getCrs( const std::string &line );
+  };
+
+} // namespace MDAL
+#endif //MDAL_MIKE21_HPP
diff --git a/external/mdal/mdal.cpp b/external/mdal/mdal.cpp
index 387b58f0a0d..6ffb6e3f8d7 100644
--- a/external/mdal/mdal.cpp
+++ b/external/mdal/mdal.cpp
@@ -21,7 +21,7 @@ static const char *EMPTY_STR = "";
 
 const char *MDAL_Version()
 {
-  return "1.0.3";
+  return "1.1.0";
 }
 
 MDAL_Status MDAL_LastStatus()
diff --git a/external/mdal/mdal_driver_manager.cpp b/external/mdal/mdal_driver_manager.cpp
index 2d99935c81f..c9c4d5ca967 100644
--- a/external/mdal/mdal_driver_manager.cpp
+++ b/external/mdal/mdal_driver_manager.cpp
@@ -11,6 +11,7 @@
 #include "frmts/mdal_binary_dat.hpp"
 #include "frmts/mdal_selafin.hpp"
 #include "frmts/mdal_esri_tin.hpp"
+#include "frmts/mdal_mike21.hpp"
 #include "frmts/mdal_dynamic_driver.hpp"
 #include "mdal_utils.hpp"
 
@@ -262,6 +263,8 @@ MDAL::DriverManager::DriverManager()
   mDrivers.push_back( std::make_shared<MDAL::DriverXdmf>() );
 #endif
 
+  mDrivers.push_back( std::make_shared<MDAL::DriverMike21>() );
+
   loadDynamicDrivers();
 }
 
diff --git a/src/providers/mdal/CMakeLists.txt b/src/providers/mdal/CMakeLists.txt
index ad87a66d8c3..10e32fdd5d4 100644
--- a/src/providers/mdal/CMakeLists.txt
+++ b/src/providers/mdal/CMakeLists.txt
@@ -48,6 +48,7 @@ if (WITH_INTERNAL_MDAL)
     ${CMAKE_SOURCE_DIR}/external/mdal/frmts/mdal_esri_tin.cpp
     ${CMAKE_SOURCE_DIR}/external/mdal/frmts/mdal_xms_tin.cpp
     ${CMAKE_SOURCE_DIR}/external/mdal/frmts/mdal_ply.cpp
+    ${CMAKE_SOURCE_DIR}/external/mdal/frmts/mdal_mike21.cpp
     ${CMAKE_SOURCE_DIR}/external/mdal/3rdparty/libplyxx/libplyxx.cpp
   )
 
@@ -67,6 +68,7 @@ if (WITH_INTERNAL_MDAL)
     ${CMAKE_SOURCE_DIR}/external/mdal/frmts/mdal_esri_tin.hpp
     ${CMAKE_SOURCE_DIR}/external/mdal/frmts/mdal_xms_tin.hpp
     ${CMAKE_SOURCE_DIR}/external/mdal/frmts/mdal_ply.hpp
+    ${CMAKE_SOURCE_DIR}/external/mdal/frmts/mdal_mike21.hpp
     ${CMAKE_SOURCE_DIR}/external/mdal/3rdparty/libplyxx.h
     ${CMAKE_SOURCE_DIR}/external/mdal/3rdparty/libplyxx/libplyxx_internal.h
     ${CMAKE_SOURCE_DIR}/external/mdal/3rdparty/textio.h