Вычисление факториала на Си

Как написать программу вычисления факториала в Си? Точно также как на любом другом языке программирования. Если число n! небольшое, то можно использовать и длинное целое. А если n=50, 100, 500 и более? То надо использовать умножение длинных чисел. Немного ее модифицировав. Так как один из множителей входит в диапазон целых чисел, то для хранения большого второго множителя используем массив целых чисел. Какой должен быть при этом размер масива? Можно оценить, используя обычный калькулятор или написав программу с дробными числами. Показатель степени в экспоненциальной форме даст граничное значение. Так для 500! потребуется 1137 цифр, а для 1000! - уже 2567. Возьмем массив до 10000 элементов. Нулевой элемент хранит последнюю цифру факториала. Тогда a[0]=1. В нашем случае нужно ораничить массив до первой цифры. Пусть это будет k. И в начале оно равно 0. По мере того, как число увеличивается, k тоже будет увеличиваться. Так для n=7 получим следующие шаги цикла по i:
i=0 a[0]=1, k=0
i=1 a[0]=1, k=0 - 0!=1!=1, поэтому цикл можно начать с двух.
i=2 a[0]=1*2=2, k=0
i=3 a[0]=2*3=6, k=0
i=4 a[0]=6*4=24, k=0
Т.е., если получившееся число больше 9, то нужно разбить его на разряды: a[0]=4, a[1]=2, и теперь k=1.
i=5 a[0]=5*4=20, здесь запомним цифру 2 для переноса в следующий разряд. Тогда a[0]=0, d=2,
      a[1]=2*5+2=12. Здесь опять запомним цифру 1 для переноса в следующий разряд: a[1]=2, d=1. Так как мы дошли до конечного k, то, также как и в предыдущем случае для i=4, увеличиваем его и можно еще раз пройти по тому же циклу j:
       a[2]=0*5+1=1, d=0.
i=6  a[0]=a[0]*6+d=0*6+0=0, d=0
       a[1]=a[1]*6+d=2*6+0=12, d=1, a[1]=2
       a[2]=a[2]*6+d=1*6+1=7, d=0.
i=7  a[0]=a[0]*7+d=0*7+0=0, d=0
       a[1]=a[1]*7+d=2*7+0=14, d=1, a[1]=4
       a[2]=a[2]*7+d=7*7+1=50, d=5, a[2]=0, k=k+1=3
       a[3]=a[3]*7+d=0*7+5=5, d=0
Получим во внешнем цикле по i от 2 до n с шагом 1 еще один цикл для умножения каждой цифры на число i. Пусть это будет цикл по j от 0 до k с шагом 1. Запишем тело цикла, учитывая все выше изложенные шаги и замечания:
      a[j]=a[j]*i+d; если a[j]<10, то d=0, иначе d=целое при делении a[j] на 10, a[j]=остаток деления a[j] на 10 (здесь важен именно такой порядок, так как если мы сразу переприсвоим элемент массива, то d может стать равным 0);  если j конечное, то увеличим k на единицу и еще раз пройдем цикл. Если a[j] будет больше 100, то этот блок повторится до тех пор, пока все разряды не займут свое место.
      Все - число n! готово. Осталось вывести. Так как k хранит все разряды результирующего числа, то в цикле по i от k до 0 выведем все цифры без пробела, начиная с первой цифры.

Получим следующую программу:

#include <stdio.h>
int main()
{
 int n;
 int a[10000];
 int d,i,j,m,k;
 scanf("%d",&n);
 for(i=0;i<=10000;i++) a[i]=0; //обнуляем все элементы массива, кроме первого  a[0]=1;

 k=0; //количество цифр в результирующем массиве цифр
  d=0; //остаток от деления на 10 для переноса в следующий разряд for(i=2;i<=n;i++)
 {
    for (j=0;j<=k;j++)
    {
      a[j]=a[j]*i+d; /*умножаем каждую цифру в массиве, начиная с конца (в нашем случае нумерация последней цифры в числе начинается с нуля.*/      if (a[j]<10) d=0;
      else
/*если результат больше 9, то запоминаем цифру для переноса в следующий разряд*/        {
            d=a[j]/10;
            a[j]=a[j]%10;
            if (j==k) ++k; /*если это последняя цифра, то увеличиваем количество цифр в массиве. Вот эта конструкция на Паскале не сработает, так как начальное и конечное значение в цикле for не изменяется в теле цикла.*/
        }
     }
  }
//вывод производим с начальной цифры, т.е. с к-го номера до нулевогоfor(i=k;i>=0;i--)printf("%d",a[i]);
return 0;
}

А теперь вопрос: почему при выполнении этой программы при любом n получаем 1? Где ошибка в программе? Замечу, что в самом алгоритме ошибок нет.