function mdl = model(ne, etype)

# Input  ne        Number of elements along one direction
#        etype     Element name
# Output mdl       Structure with model definition
#
# The function generates the model of the fan blade.
#
# Node and element numbers increase by one in circumferential
# direction and by row increments in longitudinal direction.
#
# ---------------------------------------------------------------

# Data (N, mm)

  L    =    304.8;  % Length
  R    =    609.6;  % Radius
  t    =    3.048;  % Thickness
  phiE =      0.5;  % Angle in radians
  E    =   206010;  % Young's modulus
  ny   =      0.3;  % Poisson's ratio
  rho  = 7.876e-9;  % Mass density

# Nodal point data

  np = ne + 1;  % Nodes in one row

# Model type and subtype

  mdl = struct("type", "solid", "subtype", "3d");

# Geometry and material

  geom = struct("t", t);
  mat  = struct("type", "iso", "E", E, "ny", ny, "rho", rho);

# Nodes

  phi = linspace(0, phiE, np);
  x   = linspace(0, L, np);
  y   = R * sin(phi);
  z   = R * cos(phi);

  k = 1;

  for n = 1 : np
      for m = 1 : np
          nodes(k) = struct("id",   k++,
                            "coor", [x(n), y(m), z(m)]);
      end
  end

  mdl.nodes = nodes;

# Elements

  k = 1;

  switch etype

   case "s4"

     enode1 = [1, np + 1, np + 2, 2];

     for n = 1 : ne
         enodes = enode1;
         for m = 1  :ne
             elem(k) = struct("id", k++, "type", etype,
                              "nodes", enodes++,
                              "geom", geom, "mat", mat);
         end
         enode1 += np;
     end

   case "s3"

     enode1 = [1, np + 1, np + 2]; enode2 = [np + 2, 2, 1];

     for n = 1 : ne
         enodea = enode1; enodeb = enode2;
         for m = 1 : ne
             elem(k) = struct("id", k++, "type", etype,
                              "nodes", enodea++,
                              "geom", geom, "mat", mat);
             elem(k) = struct("id", k++, "type", etype,
                              "nodes", enodeb++,
                              "geom", geom, "mat", mat);
         end
         enode1 += np; enode2 += np;
     end

  end

  mdl.elements = elem;

# Constraints

  for k = 1 : np
      fix(k) = struct("id", k, "dofs", 1 : 6);
  end

  mdl.constraints.prescribed = fix;

end
