function mfs_exp_mesh(fid, cmp)

# usage: mfs_exp_mesh(fid, cmp)
#
# Input  fid     File handle of output file
#        cmp     Structure with component
#
# This function exports the mesh in Gmsh MSH file format version 4.1.
#
# ------------------------------------------------------------------------

# Check arguments

  if (nargin != 2)
     print_usage();
  end

# Physical Names
# physicalTag equals the lifting surface index

  nofls = cmp.ls.nofls;
  fprintf(fid, "$PhysicalNames\n")
  fprintf(fid, "%d\n", nofls);
  for l = 1 : nofls
      fprintf(fid, "2 %d \"%s\"\n", l,
              sprintf("Lifting_Surface_%d", cmp.ls.ids(l)));
  endfor
  fprintf(fid, "$EndPhysicalNames\n");

# Entities: Each lifting surface defines one entity
# entityTag equals physicalTag equals lifting surface index

  fprintf(fid, "$Entities\n");
  fprintf(fid, "0 0 %d 0\n", nofls);

  for l = 1 : nofls
      P    = cmp.ls.surfs(l).P;
      cmin = min(P);
      P(:, 1) += cmp.ls.surfs(l).c';
      cmax = max(P);
      fprintf(fid, "%d %f %f %f %f %f %f 1 %d 0\n",
              l, cmin, cmax, l);
  endfor

  fprintf(fid, "$EndEntities\n");

# Find entities nodes and elements belong to

  nofnod = cmp.nodes.nofnod;
  nofpan = cmp.panels.nofpan;

  nodinf       = zeros(nofnod, 2);
  nodinf(:, 1) = 1 : nofnod;
  paninf       = zeros(nofpan, 2);
  paninf(:, 1) = 1 : nofpan;

  for l = 1 : nofls
      p1 = cmp.ls.surfs(l).p1;
      pe = cmp.ls.surfs(l).pend;
      paninf(p1 : pe, 2) = l;
      ixn = unique([cmp.panels.corner(:, p1 : pe)]);
      nodinf(ixn, 2) = l;
  endfor

  nodinf = sortrows(nodinf, 2);
  nodinf(:, 1) = cmp.nodes.ids(nodinf(:, 1));

# Build structure nodx

  nbl       = 0;   % counts blocks
  nnbl      = 1;   % counts nodes in block
  entityTag = 1;
  ix1       = 1;   % Position of first node in block 

  nodinf = [nodinf; [-1, -1]];  % End marker
  blocks = struct();

  for n = 2 : nofnod + 1

      if (nodinf(n, 2) == entityTag)
          nnbl++;
      else
         blocks(++nbl) = struct("entityDim", 2,
                                "entityTag", entityTag++,
                                "numNodesInBlock", nnbl,
                                "nodeTag", nodinf(ix1 : n - 1, 1));
         nnbl = 1; ix1 = n;
      endif

  endfor

  for l = 1 : nbl
      nodix = blocks(l).nodeTag;
      coor  = cmp.nodes.coor(nodix, :);
      blocks(l).xyz = coor';
      blocks(l).nodeTag = cmp.nodes.ids(nodix);
  endfor

  nodx = struct("numEntityBlocks", nbl, "numNodes", nofnod,
                "minNodeTag", cmp.nodes.ids(1),
                "maxNodeTag", cmp.nodes.ids(end),
                "blocks", blocks);

# Output nodal points

  mfs_exp_nodes(fid, nodx);

# Build structure eltx

  nelnod      = 4;
  ixnod       = [1, 2, 4, 3];
  elementType = 3;

  blocks = struct();

  for l = 1 : nofls

      p1    = cmp.ls.surfs(l).p1;
      pe    = cmp.ls.surfs(l).pend;
      nodel = cmp.panels.corner(ixnod, p1 : pe);

      if (columns(nodel) > 1)
         nodel = cmp.nodes.ids(nodel);
      else
         nodel = cmp.nodes.ids(nodel)';
      endif

      blocks(l) = struct("entityDim", 2, 
                         "entityTag", l,
                         "elementType", elementType,
                         "numElementsInBlock", pe - p1 + 1,
                         "elementTag", p1 : pe,
                         "nodeTag", nodel);
      
  endfor

  eltx = struct("numEntityBlocks", nofls,
                "numElements", nofpan,
                "minElementTag", 1, "maxElementTag", nofpan,
                "blocks", blocks);

# Output elements

  mfs_exp_elts(fid, eltx);

endfunction
