1 #if !defined HAVE_PARTITION_2FALL_DESC_H__
2 #define      HAVE_PARTITION_2FALL_DESC_H__
3 // This file is part of the FXT library.
4 // Copyright (C) 2012, 2013, 2014, 2019 Joerg Arndt
5 // License: GNU General Public License version 3 or later,
6 // see the file COPYING.txt in the main directory.
7 
8 #include "comb/is-partition-desc.h"
9 
10 #include "comb/comb-print.h"
11 #include "comb/print-partition-aa.h"
12 
13 #include "fxttypes.h"
14 
15 #include "fxtio.h"
16 #include "jjassert.h"
17 
18 
19 class partition_2fall_desc
20 // Partitions of n is a partition a[1] + a[2] + ... + a[m] = n
21 //   such that 2*a[k] <= a[k-1].
22 // Representation as weakly descending list of parts.
23 // Lexicographic order.
24 // Cf. OEIS sequence A000929.
25 // Equinumerous to s-partitions (partitions into Mersenne numbers),
26 //  cf. class partition_s_desc.
27 {
28 public:
29     ulong n_;   // integer partition of n
30     ulong m_;   // current partition has m parts
31     ulong *a_;  // partition: a[1] + a[2] + ... + a[m] = n
32 
33     partition_2fall_desc(const partition_2fall_desc&) = delete;
34     partition_2fall_desc & operator = (const partition_2fall_desc&) = delete;
35 
36 private:
mers_t(ulong s)37     ulong mers_t(ulong s)
38     // Return greatest t such that 2^t-1 <= s.
39     {
40         ulong t = 0;
41         ulong b = 1;
42         while ( s+1 > b )  { ++t;  b <<= 1; }
43 
44         if ( s < b-1 )  { --t;  b >>= 1; }
45 
46         return t;
47     }
48 
write_tail(ulong s,ulong j)49     ulong write_tail(ulong s, ulong j)
50     // Write lexicographically first partition of s, starting at index j.
51     // Return last index written to.
52     // Undefined for s == 0.
53     {
54         // First part in first partition of n is A046699(n+1)
55 #if 0
56         // todo: optimize using A046699
57 #else
58         ulong h = mers_t(s);
59         ulong t = j + h - 1;  // return value
60         ulong b = 1UL << h;
61         s -= (b-1);
62         b >>= 1;
63         for (ulong i=j; b!=0;  ++i, b>>=1)  a_[i] = b;  // write
64 
65         while ( s != 0 )
66         {
67 #if 1
68             if ( s<=2 )
69             {
70                 a_[j] += s;
71                 return t;
72             }
73 #endif
74             h = mers_t(s);
75             b = 1UL << h;
76             s -= (b-1);
77             b >>= 1;
78             for (ulong i=j; b!=0;  ++i, b>>=1)  a_[i] += b;  // add
79         }
80 
81         return t;
82 #endif
83     }
84 
85 
86 public:
partition_2fall_desc(ulong n)87     explicit partition_2fall_desc(ulong n)
88     {
89         n_ = n;
90         ulong k = 1;
91         if ( n_ != 0 )  k += mers_t(n_);
92         a_ = new ulong[k];
93         first();
94     }
95 
~partition_2fall_desc()96     ~partition_2fall_desc()
97     { delete [] a_; }
98 
data()99     const ulong * data()  const  { return  a_ + 1; }
100 
num_parts()101     ulong num_parts()  const  { return m_; }
102 
first()103     void first()
104     {
105         a_[0] = 3 * n_;  // read-only
106         m_ = 0;
107         if ( n_!=0 )  m_ = write_tail(n_, 1);
108     }
109 
110 //    void last()
111 //    {
112 //        a_[0] = 3 * n_;  // read-only
113 //        if ( n_==0 )
114 //        {
115 //            m_ = 0;
116 //        }
117 //        else
118 //        {
119 //            m_ = 1;
120 //            a_[1] = n_;
121 //        }
122 //    }
123 
124 
next()125     ulong next()
126     // Return number of parts of generated partition.
127     // Return zero if the current is the last partition.
128     {
129         if ( m_ <= 1 )  return 0;  // current is last
130 
131         ulong j = m_;
132         ulong s = a_[j];
133 
134         ulong y = a_[j-1];
135         ulong x = a_[j-2];  // can read sentinel a[0]
136 
137         --j;
138         while ( 2*(y+1) > x  )  // search for y we can inrease
139         {
140             s += y;
141             y = x;
142             --j;
143             x = a_[j-1];  // can read sentinel a[0]
144         }
145 
146         a_[j] = y + 1;
147         s -= 1;
148 
149         if ( s == 0 )
150         {
151             --m_;
152             return m_;
153         }
154 
155         m_ = write_tail(s, j+1);
156         return m_;
157     }
158 
159 //    ulong prev()
160 //    // Return number of parts of generated partition.
161 //    // Return zero if the current is the first partition.
162 //    {
163 //        ulong z = a_[m_];
164 //        if ( z >= 3 )  // split off one unit
165 //        {
166 //            cout << " :: A 1" << endl;
167 //            a_[m_] = z - 1;
168 //            ++m_;
169 //            a_[m_] = 1;
170 //            return m_;
171 //        }
172 //
173 //        if ( n_ <= 2 )  return 0;
174 //
175 //        ulong j = m_ - 1;
176 //        ulong y = a_[j];
177 //        ulong s = z;
178 //
179 //        if ( 2*(s + 1) <= y - 1 )
180 //        {
181 //            cout << " :: A 2" << endl;
182 //            a_[m_-1] = y - 1;
183 //            a_[m_] = s + 1;
184 //            return m_;
185 //        }
186 //
187 //        // OK up to here
188 //
189 //        return 0;
190 //    }
191 
192     void print(const char *bla, bool dfz=false)  const
193     { print_vec(bla, data(), num_parts(), dfz); }
194 
print_aa()195     void print_aa()  const  // ASCII art
196     { print_partition_desc_aa(data(), m_); }
197 
print_conj_aa()198     void print_conj_aa()  const  // ASCII art
199     { print_partition_desc_conj_aa(data(), m_); }
200 
201 
OK()202     bool OK()  const
203     {
204         if ( ! is_partition_desc(data(), num_parts(), n_) )  return false;
205 
206         for (ulong j=2; j<=m_; ++j)  // 2*a[j] <= a[j-1] ?
207             if ( 2*a_[j] > a_[j-1] )  return false;
208 
209         return true;
210     }
211 };
212 // -------------------------
213 
214 
215 #endif  // !defined HAVE_PARTITION_2FALL_DESC_H__
216