function [nxu,nyu,nzu,nxv,nyv,nzv,nxw,nyw,nzw,...
          xu1,yu1,zu1,xv1,yv1,zv1,xw1,yw1,zw1,xp1,yp1,zp1,dx,dy,dz]=...
            generate_grid(a,b,c,d,e,f,nxp,nyp,nzp,x_periodic,y_periodic,z_periodic)
  
  if x_periodic==1
    nxu=nxp;
    nxv=nxp;
    nxw=nxp;
  else
    nxu=nxp+1;
    nxv=nxp;
    nxw=nxp;
  end
  if y_periodic==1
    nyu=nyp;
    nyv=nyp;
    nyw=nyp;
  else
    nyu=nyp;
    nyv=nyp+1;
    nyw=nyp;
  end
  if z_periodic==1
    nzu=nzp;
    nzv=nzp;
    nzw=nzp;
  else
    nzu=nzp;
    nzv=nzp;
    nzw=nzp+1;
  end
  %grid step:
  if x_periodic==1
    dx=(b-a)/(nxp);
  else
    dx=(b-a)/(nxp-1);
  end
  if y_periodic==1
    dy=(d-c)/(nyp);
  else
    dy=(d-c)/(nyp-1);
  end
  if z_periodic==1
    dz=(f-e)/(nzp);
  else
    dz=(f-e)/(nzp-1);
  end
  %
  if x_periodic==1
    xu1=a+((1:nxu)-1)*dx;
  else
    xu1=a+((1:nxu)-3/2)*dx;
  end
  if y_periodic==1
    yu1=c+((1:nyu)-1/2)*dy;
  else
    yu1=c+((1:nyu)-1)*dy;
  end
  if z_periodic==1
    zu1=e+((1:nzu)-1/2)*dz;
  else
    zu1=e+((1:nzu)-1)*dz;
  end
  %
  if x_periodic==1
    xv1=a+((1:nxv)-1/2)*dx;
  else
    xv1=a+((1:nxv)-1)*dx;
  end
  if y_periodic==1
    yv1=c+((1:nyv)-1)*dy;
  else
    yv1=c+((1:nyv)-3/2)*dy;
  end
  if z_periodic==1
    zv1=e+((1:nzv)-1/2)*dz;
  else
    zv1=e+((1:nzv)-1)*dz;
  end
  %
  if x_periodic==1
    xw1=a+((1:nxw)-1/2)*dx;
  else
    xw1=a+((1:nxw)-1)*dx;
  end
  if y_periodic==1
    yw1=c+((1:nyw)-1/2)*dy;
  else
    yw1=c+((1:nyw)-1)*dy;
  end
  if z_periodic==1
    zw1=e+((1:nzw)-1)*dz;
  else
    zw1=e+((1:nzw)-3/2)*dz;
  end
  %
  if x_periodic==1
    xp1=a+((1:nxp)-1/2)*dx;
  else
    xp1=a+((1:nxp)-1)*dx;
  end
  if y_periodic==1
    yp1=c+((1:nyp)-1/2)*dy;
  else
    yp1=c+((1:nyp)-1)*dy;
  end
  if z_periodic==1
    zp1=e+((1:nzp)-1/2)*dz;
  else
    zp1=e+((1:nzp)-1)*dz;
  end
end