← Back to branch summary

~ubuntu-branches/debian/sid/pgadmin3/sid

« back to all changes in this revision

Viewing changes to pgadmin/pgscript/utilities/m_apm/mapmfact.cpp

  • Committer: Bazaar Package Importer
  • Author(s): Gerfried Fuchs
  • Date: 2009-07-30 12:27:16 UTC
  • mfrom: (1.1.6 upstream)
  • Revision ID: james.westby@ubuntu.com-20090730122716-fddbh42on721bbs2
Tags: 1.10.0-1
* New upstream release.
* Adjusted watch file to match release candidates.
* Updated to Standards-Version 3.8.2:
  - Moved to Section: database.
  - Add DEB_BUILD_OPTIONS support for parallel building.
  - Move from findstring to filter suggestion for DEB_BUILD_OPTIONS parsing.
* pgagent got split into its own separate source package by upstream.
* Exclude Docs.vcproj from installation.
* Move doc-base.enus from pgadmin3 to pgadmin3-data package, the files are
  in there too.

Show diffs side-by-side

added added

removed removed

Lines of Context:
pgadmin/pgscript/utilities/m_apm/mapmfact.cpp
 
1
 
 
2
/* 
 
3
 *  M_APM  -  mapmfact.c
 
4
 *
 
5
 *  Copyright (C) 1999 - 2007   Michael C. Ring
 
6
 *
 
7
 *  Permission to use, copy, and distribute this software and its
 
8
 *  documentation for any purpose with or without fee is hereby granted,
 
9
 *  provided that the above copyright notice appear in all copies and
 
10
 *  that both that copyright notice and this permission notice appear
 
11
 *  in supporting documentation.
 
12
 *
 
13
 *  Permission to modify the software is granted. Permission to distribute
 
14
 *  the modified code is granted. Modifications are to be distributed by
 
15
 *  using the file 'license.txt' as a template to modify the file header.
 
16
 *  'license.txt' is available in the official MAPM distribution.
 
17
 *
 
18
 *  This software is provided "as is" without express or implied warranty.
 
19
 */
 
20
 
 
21
/*
 
22
 *
 
23
 *      This file contains the FACTORIAL function.
 
24
 *
 
25
 */
 
26
 
 
27
/*
 
28
 *      Brief explanation of the factorial algorithm.
 
29
 *      ----------------------------------------------
 
30
 *
 
31
 *      The old algorithm simply multiplied N * (N-1) * (N-2) etc, until 
 
32
 *      the number counted down to '2'. So one term of the multiplication 
 
33
 *      kept getting bigger while multiplying by the next number in the 
 
34
 *      sequence. 
 
35
 *
 
36
 *      The new algorithm takes advantage of the fast multiplication 
 
37
 *      algorithm. The "ideal" setup for fast multiplication is when 
 
38
 *      both numbers have approx the same number of significant digits 
 
39
 *      and the number of digits is very near (but not over) an exact 
 
40
 *      power of 2.
 
41
 *
 
42
 *      So, we will multiply N * (N-1) * (N-2), etc until the number of
 
43
 *      significant digits is approx 256.
 
44
 *
 
45
 *      Store this temp product into an array.
 
46
 *
 
47
 *      Then we will multiply the next sequence until the number of
 
48
 *      significant digits is approx 256.
 
49
 *
 
50
 *      Store this temp product into the next element of the array.
 
51
 *
 
52
 *      Continue until we've counted down to 2.
 
53
 *
 
54
 *      We now have an array of numbers with approx the same number
 
55
 *      of digits (except for the last element, depending on where it 
 
56
 *      ended.) Now multiply each of the array elements together to
 
57
 *      get the final product.
 
58
 *
 
59
 *      The array multiplies are done as follows (assume we used 11
 
60
 *      array elements for this example, indicated by [0] - [10] ) :
 
61
 *
 
62
 *      initial    iter-1     iter-2       iter-3     iter-4
 
63
 *
 
64
 *        [0] 
 
65
 *           *  ->  [0]
 
66
 *        [1]
 
67
 *                      * ->    [0]
 
68
 *
 
69
 *        [2] 
 
70
 *           *  ->  [1]
 
71
 *        [3]
 
72
 *                                   * ->   [0] 
 
73
 *
 
74
 *        [4] 
 
75
 *           *  ->  [2]
 
76
 *        [5]
 
77
 *
 
78
 *                      * ->    [1]
 
79
 *
 
80
 *        [6] 
 
81
 *           *  ->  [3]                           *  ->  [0]
 
82
 *        [7]
 
83
 *
 
84
 *
 
85
 *        [8] 
 
86
 *           *  ->  [4]
 
87
 *        [9]
 
88
 *                      * ->    [2]    ->   [1]
 
89
 *
 
90
 *
 
91
 *        [10]  ->  [5]
 
92
 *
 
93
 */
 
94
 
 
95
#include "pgAdmin3.h"
 
96
#include "pgscript/utilities/mapm-lib/m_apm_lc.h"
 
97
 
 
98
/* define size of local array for temp storage */
 
99
 
 
100
#define NDIM 32
 
101
 
 
102
/****************************************************************************/
 
103
void    m_apm_factorial(M_APM moutput, M_APM minput)
 
104
{
 
105
int     ii, nmul, ndigits, nd, jj, kk, mm, ct;
 
106
M_APM   array[NDIM];
 
107
M_APM   iprod1, iprod2, tmp1, tmp2;
 
108
 
 
109
/* return 1 for any input <= 1 */
 
110
 
 
111
if (m_apm_compare(minput, MM_One) <= 0)
 
112
  {
 
113
   m_apm_copy(moutput, MM_One);
 
114
   return;
 
115
  }
 
116
 
 
117
ct       = 0;
 
118
mm       = NDIM - 2;
 
119
ndigits  = 256;
 
120
nd       = ndigits - 20;
 
121
tmp1     = m_apm_init();
 
122
tmp2     = m_apm_init();
 
123
iprod1   = m_apm_init();
 
124
iprod2   = m_apm_init();
 
125
array[0] = m_apm_init();
 
126
 
 
127
m_apm_copy(tmp2, minput);
 
128
 
 
129
/* loop until multiply count-down has reached '2' */
 
130
 
 
131
while (TRUE)
 
132
  {
 
133
   m_apm_copy(iprod1, MM_One);
 
134
 
 
135
   /* 
 
136
    *   loop until the number of significant digits in this 
 
137
    *   partial result is slightly less than 256
 
138
    */
 
139
 
 
140
   while (TRUE)
 
141
     {
 
142
      m_apm_multiply(iprod2, iprod1, tmp2);
 
143
 
 
144
      m_apm_subtract(tmp1, tmp2, MM_One);
 
145
 
 
146
      m_apm_multiply(iprod1, iprod2, tmp1);
 
147
 
 
148
      /*
 
149
       *  I know, I know.  There just isn't a *clean* way 
 
150
       *  to break out of 2 nested loops.
 
151
       */
 
152
 
 
153
      if (m_apm_compare(tmp1, MM_Two) <= 0)
 
154
        goto PHASE2;
 
155
 
 
156
      m_apm_subtract(tmp2, tmp1, MM_One);
 
157
 
 
158
      if (iprod1->m_apm_datalength > nd)
 
159
        break;
 
160
     }
 
161
 
 
162
   if (ct == (NDIM - 1))
 
163
     {
 
164
      /* 
 
165
       *    if the array has filled up, start multiplying
 
166
       *    some of the partial products now.
 
167
       */
 
168
 
 
169
      m_apm_copy(tmp1, array[mm]);
 
170
      m_apm_multiply(array[mm], iprod1, tmp1);
 
171
 
 
172
      if (mm == 0)
 
173
        {
 
174
         mm = NDIM - 2;
 
175
         ndigits = ndigits << 1;
 
176
         nd = ndigits - 20;
 
177
        }
 
178
      else
 
179
         mm--;
 
180
     }
 
181
   else
 
182
     {
 
183
      /* 
 
184
       *    store this partial product in the array
 
185
       *    and allocate the next array element
 
186
       */
 
187
 
 
188
      m_apm_copy(array[ct], iprod1);
 
189
      array[++ct] = m_apm_init();
 
190
     }
 
191
  }
 
192
 
 
193
PHASE2:
 
194
 
 
195
m_apm_copy(array[ct], iprod1);
 
196
 
 
197
kk = ct;
 
198
 
 
199
while (kk != 0)
 
200
  {
 
201
   ii = 0;
 
202
   jj = 0;
 
203
   nmul = (kk + 1) >> 1;
 
204
 
 
205
   while (TRUE)
 
206
     {
 
207
      /* must use tmp var when ii,jj point to same element */
 
208
 
 
209
      if (ii == 0)
 
210
        {
 
211
         m_apm_copy(tmp1, array[ii]);
 
212
         m_apm_multiply(array[jj], tmp1, array[ii+1]);
 
213
        }
 
214
      else
 
215
         m_apm_multiply(array[jj], array[ii], array[ii+1]);
 
216
 
 
217
      if (++jj == nmul)
 
218
        break;
 
219
 
 
220
      ii += 2;
 
221
     }
 
222
 
 
223
   if ((kk & 1) == 0)
 
224
     {
 
225
      jj = kk >> 1;
 
226
      m_apm_copy(array[jj], array[kk]);
 
227
     }
 
228
 
 
229
   kk = kk >> 1;
 
230
  }
 
231
 
 
232
m_apm_copy(moutput, array[0]);
 
233
 
 
234
for (ii=0; ii <= ct; ii++)
 
235
  {
 
236
   m_apm_free(array[ii]);
 
237
  }
 
238
 
 
239
m_apm_free(tmp1);
 
240
m_apm_free(tmp2);
 
241
m_apm_free(iprod1);
 
242
m_apm_free(iprod2);
 
243
}
 
244
/****************************************************************************/