Класс больших чисел. Деление по Кнуту.

@Ultrator · Регистрация: 28.04.2009

Студворк — интернет-сервис помощи студентам

Кнут т2 стр. 300 алгоритм D. Кто реализовал для b=2^32 или 2^16 - поделитесь, плз, если не жалко !
(А то у меня много неясных моментов. Да уже и сделал - но работает неправильно, я не стал разбираться, стёр всё нафиг.) Плз, хелп!

06.03.2010, 12:52

За основу брал реализацию из книги "Алгоритмические трюки для программиста"

C++
1
2
3
4
5
6
7
8
9
10
11
int nlz(unsigned int x)
{
    unsigned int y;
    int n=32;
    y=x>>16;if(y){n-=16;x=y;}
    y=x>>8; if(y){n-=8;x=y;}
    y=x>>4; if(y){n-=4;x=y;}
    y=x>>2; if(y){n-=2;x=y;}
    y=x>>1; if(y)return n-2;
    return n-x;
}

Основание 2^16, прямой порядок следования:

C++
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
int divmnu(unsigned short q[],unsigned short r[],const unsigned short u[],const unsigned short v[],int m,int n)
{
    const unsigned int b=0x10000; // Основание чисел (16 битов)
    unsigned short *un,*vn; // Нормализованный вид u и v
    unsigned int qhat; // Предполагаемая цифра частного
    unsigned int rhat; // Остаток
    unsigned int p; // Произведение двух цифр
    int s,i,j,t,k;
    if(m<n || n<=0 || v[n-1]==0)
        return 1; // Возвращается при некорректном параметре
    if(n==1) // Частный случай делителя из одной цифры
    {
        k=0;
        for(j=m-1;j>=0;j--)
        {
            q[j]=(k*b + u[j])/v[0];
            k=(k*b + u[j])-q[j]*v[0];
        }
        if(r!=0)r[0]=k;
        return 0;
    }
// Нормализация путем сдвига v влево, такого, что
// старший бит становится единичным, и сдвига u влево
// на ту же величину. Нам может потребоваться добавить
// старшую цифру к частному; мы делаем это безусловно
    s=nlz(v[n-1])-16; // 0 <= s <= 16.
    vn=(unsigned short *)alloca(2*n);
    for(i=n-1;i>0;i--)
        vn[i]=(v[i]<<s) | (v[i-1]>>16-s);
    vn[0]=v[0]<<s;
    un=(unsigned short *)alloca(2*(m+1));
    un[m]=u[m-1]>>16-s;
    for(i=m-1;i>0;i--)
        un[i]=(u[i]<<s) | (u[i-1]>>16-s);
    un[0]=u[0]<<s;
    for(j=m-n;j>=0;j--) // Главный цикл
    { // Вычисляем оценку q[j]
        qhat=(un[j+n]*b + un[j+n-1])/vn[n-1];
        rhat=(un[j+n]*b + un[j+n-1])-qhat*vn[n-1];
again:
        if(qhat>=b || qhat*vn[n-2]>b*rhat + un[j+n-2])
        {
            qhat=qhat-1;
            rhat=rhat+vn[n-1];
            if(rhat<b)goto again;
        }
// Умножение и вычитание
        k=0;
        for(i=0;i<n;i++)
        {
            p=qhat*vn[i];
            t=un[i+j]-k-(p&0xFFFF);
            un[i+j]=t;
            k=(p>>16)-(t>>16);
        }
        t=un[j+n]-k;
        un[j+n]=t;
        q[j]=qhat; // Сохранение цифры частного
        if(t<0) // Если мы вычли слишком много,
        { // вернем назад
            q[j]=q[j]-1;
            k=0;
            for(i=0;i<n;i++)
            {
                t=un[i+j]+vn[i]+k;
                un[i+j]=t;
                k=t>>16;
            }
            un[j+n]=un[j+n]+k;
        }
    } // for j
// Если вызывающей функции нужно значение остатка,
// денормализуем и возвращаем его
    if(r!=0)
    {
        for(i=0;i<n;i++)
        r[i]=(un[i]>>s) | (un[i+1]<<16-s);
    }
    return 0;
}

Основание 2^32, прямой порядок следования:

C++
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
int divmnu(unsigned int q[],unsigned int r[],const unsigned int u[],const unsigned int v[],int m,int n)
{
    const unsigned _int64 b=0x100000000; // Основание чисел (32 бита)
    unsigned int *un,*vn; // Нормализованный вид u и v
    int i,j;
    unsigned int k;
    if(m<n || n<=0 || v[n-1]==0)
        return 1; // Возвращается при некорректном параметре
    if(n==1) // Частный случай делителя из одной цифры
    {
        k=0;
        for(j=m-1;j>=0;j--)
        {
            q[j]=(k*b + u[j])/v[0];
            k=(k*b + u[j])-q[j]*v[0];
        }
        if(r!=0)r[0]=k;
        return 0;
    }
// Нормализация путем сдвига v влево, такого, что
// старший бит становится единичным, и сдвига u влево
// на ту же величину. Нам может потребоваться добавить
// старшую цифру к частному; мы делаем это безусловно
    int s=nlz(v[n-1]); // 0 <= s <= 32.
    vn=(unsigned int *)alloca(4*n);
    for(i=n-1;i>0;i--)
        vn[i]=(v[i]<<s) | (v[i-1]>>32-s);
    vn[0]=v[0]<<s;
    un=(unsigned int *)alloca(4*(m+1));
    un[m]=u[m-1]>>32-s;
    for(i=m-1;i>0;i--)
        un[i]=(u[i]<<s) | (u[i-1]>>32-s);
    un[0]=u[0]<<s;
    for(j=m-n;j>=0;j--) // Главный цикл
    { // Вычисляем оценку q[j]
        unsigned int qhat=(un[j+n]*b + un[j+n-1])/vn[n-1]; // Предполагаемая цифра частного
        unsigned _int64 rhat=(un[j+n]*b + un[j+n-1])-qhat*vn[n-1]; // Остаток
        while(1)
        {
            if(qhat>=b || qhat*vn[n-2]>b*rhat + un[j+n-2])
            {
                qhat=qhat-1;
                rhat=rhat+vn[n-1];
                if(rhat<b)continue;
            }
            break;
        }
// Умножение и вычитание
        k=0;
        unsigned _int64 t;
        for(i=0;i<n;i++)
        {
            unsigned _int64 p=(unsigned _int64)qhat*vn[i];
            t=un[i+j]-k-(p&0xFFFFFFFF);
            un[i+j]=t;
            k=(p>>32)-(t>>32);
        }
        t=un[j+n]-k;
        un[j+n]=t;
        q[j]=qhat; // Сохранение цифры частного
        if(t<0) // Если мы вычли слишком много,
        { // вернем назад
            q[j]=q[j]-1;
            k=0;
            for(i=0;i<n;i++)
            {
                t=un[i+j]+vn[i]+k;
                un[i+j]=t;
                k=t>>32;
            }
            un[j+n]=un[j+n]+k;
        }
    } // for j
// Если вызывающей функции нужно значение остатка,
// денормализуем и возвращаем его
    if(r!=0)
    {
        for(i=0;i<n;i++)
        r[i]=(un[i]>>s) | (un[i+1]<<32-s);
    }
    return 0;
}

Основание 2^32, обратный порядок следования:

C++
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
int divmnu(unsigned int q[],unsigned int r[],const unsigned int u[],const unsigned int v[],int m,int n)
{
    const unsigned _int64 b=0x100000000; // Основание чисел (32 бита)
    unsigned int *un,*vn; // Нормализованный вид u и v
    int i,j;
    unsigned int k;
    if(m<n || n<=0 || v[0]==0)
        return 1; // Возвращается при некорректном параметре
    if(n==1) // Частный случай делителя из одной цифры
    {
        k=0;
        for(j=0;j<m;j++)
        {
            q[j]=(k*b + u[j])/v[n-1];
            k=(k*b + u[j])-q[j]*v[n-1];
        }
        if(r!=0)r[0]=k;
        return 0;
    }
// Нормализация путем сдвига v влево, такого, что
// старший бит становится единичным, и сдвига u влево
// на ту же величину. Нам может потребоваться добавить
// старшую цифру к частному; мы делаем это безусловно
    int s=nlz(v[0]); // 0 <= s <= 32.
    vn=(unsigned int *)alloca(4*n);
    for(i=0;i<n-1;i++)
        vn[i]=(v[i]<<s) | (v[i+1]>>32-s);
    vn[n-1]=v[n-1]<<s;
    un=(unsigned int *)alloca(4*(m+1));
    un[0]=u[0]>>32-s;
    for(i=1;i<m;i++)
        un[i]=(u[i-1]<<s) | (u[i]>>32-s);
    un[m]=u[m-1]<<s;
    for(j=0;j<=m-n;j++) // Главный цикл
    { // Вычисляем оценку q[j]
        unsigned int qhat=(un[j]*b + un[j+1])/vn[0]; // Предполагаемая цифра частного
        unsigned _int64 rhat=(un[j]*b + un[j+1])-qhat*vn[0]; // Остаток
        while(1)
        {
            if(qhat>=b || qhat*vn[1]>b*rhat + un[j+1])
            {
                qhat=qhat-1;
                rhat=rhat+vn[0];
                if(rhat<b)continue;
            }
            break;
        }
// Умножение и вычитание
        k=0;
        unsigned _int64 t;
        for(i=n;i>0;i--)
        {
            unsigned _int64 p=(unsigned _int64)qhat*vn[i-1];
            t=un[i+j]-k-(p&0xFFFFFFFF);
            un[i+j]=t;
            k=(p>>32)-(t>>32);
        }
        t=un[j]-k;
        un[j]=t;
        q[j]=qhat; // Сохранение цифры частного
        if(t<0) // Если мы вычли слишком много,
        { // вернем назад
            q[j]=q[j]-1;
            k=0;
            for(i=n;i>0;i--)
            {
                t=un[i+j]+vn[i-1]+k;
                un[i+j]=t;
                k=t>>32;
            }
            un[j]=un[j]+k;
        }
    } // for j
// Если вызывающей функции нужно значение остатка,
// денормализуем и возвращаем его
    if(r!=0)
    {
        for(i=n-1;i>=0;i--)
        r[i]=(un[i+1]>>s) | (un[i]<<32-s);
    }
    return 0;
}

Добавлено через 23 минуты
В денормализации косяк обнаружился, для обратного порядка следования исправил:

C++
1
2
for(i=m-1;i>0;i--)
r[i-1]=(un[i+1]>>s) | (un[i]<<32-s);

Добавлено через 8 минут
Не, опять не правильно, вроде так:

C++
1
for(i=n;i>0;i--)

Новые блоги и статьи Все статьи Все блоги /
PhpStorm 2025.3: WSL Terminal всегда стартует в ~ and_y87 14.12.2025 PhpStorm 2025. 3: WSL Terminal всегда стартует в ~ (home), игнорируя директорию проекта Симптом: После обновления до PhpStorm 2025. 3 встроенный терминал WSL открывается в домашней директории. . .	Access VikBal 11.12.2025 Помогите пожалуйста !! Как объединить 2 одинаковые БД Access с разными данными.	Новый ноутбук volvo 07.12.2025 Всем привет. По скидке в "черную пятницу" взял себе новый ноутбук Lenovo ThinkBook 16 G7 на Амазоне: Ryzen 5 7533HS 64 Gb DDR5 1Tb NVMe 16" Full HD Display Win11 Pro	Музыка, написанная Искусственным Интеллектом volvo 04.12.2025 Всем привет. Некоторое время назад меня заинтересовало, что уже умеет ИИ в плане написания музыки для песен, и, собственно, исполнения этих самых песен. Стихов у нас много, уже вышли 4 книги, еще 3. . .	От async/await к виртуальным потокам в Python IndentationError 23.11.2025 Армин Ронахер поставил под сомнение async/ await. Создатель Flask заявляет: цветные функции - провал, виртуальные потоки - решение. Не threading-динозавры, а новое поколение лёгких потоков. Откат?. . .
Поиск "дружественных имён" СОМ портов Argus19 22.11.2025 Поиск "дружественных имён" СОМ портов На странице: https:/ / norseev. ru/ 2018/ 01/ 04/ comportlist_windows/ нашёл схожую тему. Там приведён код на С++, который показывает только имена СОМ портов, типа,. . .	Сколько Государство потратило денег на меня, обеспечивая инсулином. Programma_Boinc 20.11.2025 Сколько Государство потратило денег на меня, обеспечивая инсулином. Вот решила сделать интересный приблизительный подсчет, сколько государство потратило на меня денег на покупку инсулинов. . . .	Ломающие изменения в C#.NStar Alpha Etyuhibosecyu 20.11.2025 Уже можно не только тестировать, но и пользоваться C#. NStar - писать оконные приложения, содержащие надписи, кнопки, текстовые поля и даже изображения, например, моя игра "Три в ряд" написана на этом. . .	Мысли в слух kumehtar 18.11.2025 Кстати, совсем недавно имел разговор на тему медитаций с людьми. И обнаружил, что они вообще не понимают что такое медитация и зачем она нужна. Самые базовые вещи. Для них это - когда просто люди. . .	Создание Single Page Application на фреймах krapotkin 16.11.2025 Статья исключительно для начинающих. Подходы оригинальностью не блещут. В век Веб все очень привыкли к дизайну Single-Page-Application . Быстренько разберем подход "на фреймах". Мы делаем одну. . .

@Ultrator 14 / 10 / 1 Регистрация: 28.04.2009 Сообщений: 219

	Класс больших чисел. Деление по Кнуту. 07.02.2010, 16:25. Показов 6324. Ответов 1 Метки нет (Все метки) Кнут т2 стр. 300 алгоритм D. Кто реализовал для b=2^32 или 2^16 - поделитесь, плз, если не жалко ! (А то у меня много неясных моментов. Да уже и сделал - но работает неправильно, я не стал разбираться, стёр всё нафиг.) Плз, хелп! 0

Класс больших чисел. Деление по Кнуту.

Решение