void red(int iz1, int iz2, int jz1, int jz2, int jm1, int jm2,
		 int jmf, int ic1, int jc1, int jcf, int kc, double c[],
		 int nci, int ncj, int nck, double s[], int nsi, int nsj)
{
    int loff = jc1 - jm1;
    int j,l,i,ic = ic1;
	double vx;
    for (j = jz1; j<=jz2; j++)
	{
        for (l = jm1; l<=jm2; l++)
		{
            vx = c[((ic - 1) * ncj + l + loff - 1) * nck + kc];
            for (i = iz1; i<=iz2; i++)
			{
               s[(i-1)*nsj+l]=s[(i-1)*nsj+l]-s[(i-1)*nsj+j]*vx;
            }
        }
        vx = c[((ic - 1) * ncj + jcf - 1) * nck + kc];
        for (i = iz1; i<=iz2; i++)
		{
            s[(i-1)*nsj+jmf]=s[(i-1)*nsj+jmf]-s[(i-1)*nsj+j]*vx;
        }
        ic = ic + 1;
    }
}