@@ -180,10 +180,18 @@ static inline void bc_standard_sqrt(bc_num *num, size_t rscale, size_t num_calc_
180
180
}
181
181
guess_vector [guess_arr_size - 1 ] = 0 ;
182
182
183
- size_t quot_size = n_arr_size - (guess_arr_size - 1 ) + 1 ;
184
-
185
183
BC_VECTOR two [1 ] = { 2 };
186
184
185
+ /* The precision (number of vectors) used for the calculation.
186
+ * Since the initial value uses two vectors, the initial precision is set to 2. */
187
+ size_t guess_precision = 2 ;
188
+ size_t guess_offset = guess_arr_size - 1 - guess_precision ;
189
+ size_t n_offset = guess_offset * 2 ;
190
+ size_t n_precision = n_arr_size - n_offset ;
191
+ size_t quot_size = n_precision - (guess_precision ) + 1 ;
192
+ size_t guess_use_len = guess_top_vector_len + BC_VECTOR_SIZE ;
193
+ bool updated_precision = false;
194
+
187
195
/**
188
196
* Newton's algorithm. Iterative expression is `x_{n+1} = (x_n + a / x_n) / 2`
189
197
* If break down the calculation into detailed steps, it looks like this:
@@ -194,16 +202,25 @@ static inline void bc_standard_sqrt(bc_num *num, size_t rscale, size_t num_calc_
194
202
*/
195
203
bool done = false;
196
204
do {
205
+ if (updated_precision ) {
206
+ guess_offset = guess_arr_size - 1 - guess_precision ;
207
+ n_offset = guess_offset * 2 ;
208
+ n_precision = n_arr_size - n_offset ;
209
+ quot_size = n_precision - (guess_precision ) + 1 ;
210
+ guess_use_len = guess_top_vector_len + (guess_precision - 1 ) * BC_VECTOR_SIZE ;
211
+ updated_precision = false;
212
+ }
213
+
197
214
/* Since the value changes during division by successive approximation, use a copied version of it. */
198
- for (size_t i = 0 ; i < n_arr_size ; i ++ ) {
215
+ for (size_t i = n_offset ; i < n_arr_size ; i ++ ) {
199
216
n_vector_copy [i ] = n_vector [i ];
200
217
}
201
218
202
219
/* 1. quot = a / x_n */
203
220
bool div_ret = bc_divide_vector (
204
- n_vector_copy , n_arr_size ,
205
- guess_vector , guess_arr_size - 1 , guess_full_len ,
206
- tmp_div_ret_vector , quot_size
221
+ n_vector_copy + n_offset , n_precision ,
222
+ guess_vector + guess_offset , guess_precision , guess_use_len ,
223
+ tmp_div_ret_vector + guess_offset , quot_size
207
224
);
208
225
ZEND_ASSERT (div_ret );
209
226
@@ -213,7 +230,7 @@ static inline void bc_standard_sqrt(bc_num *num, size_t rscale, size_t num_calc_
213
230
214
231
/* 2. add = x_n + quot1 */
215
232
int carry = 0 ;
216
- for (size_t i = 0 ; i < guess_arr_size - 1 ; i ++ ) {
233
+ for (size_t i = guess_offset ; i < guess_arr_size ; i ++ ) {
217
234
guess_vector [i ] = guess1_vector [i ] + tmp_div_ret_vector [i ] + carry ;
218
235
if (guess_vector [i ] >= BC_VECTOR_BOUNDARY_NUM ) {
219
236
guess_vector [i ] -= BC_VECTOR_BOUNDARY_NUM ;
@@ -226,27 +243,37 @@ static inline void bc_standard_sqrt(bc_num *num, size_t rscale, size_t num_calc_
226
243
227
244
/* 3. x_{n+1} = add / 2 */
228
245
div_ret = bc_divide_vector (
229
- guess_vector , guess_arr_size ,
246
+ guess_vector + guess_offset , guess_precision + 1 ,
230
247
two , 1 , 1 ,
231
- tmp_div_ret_vector , guess_arr_size
248
+ tmp_div_ret_vector + guess_offset , guess_precision + 1
232
249
);
233
250
ZEND_ASSERT (div_ret );
234
251
235
- for (size_t i = 0 ; i < guess_arr_size ; i ++ ) {
252
+ for (size_t i = guess_offset ; i < guess_arr_size ; i ++ ) {
236
253
guess_vector [i ] = tmp_div_ret_vector [i ];
237
254
}
238
255
239
256
/* 4. repeat until the difference between the `x_n` and `x_{n+1}` is less than or equal to 1. */
240
- size_t diff = guess_vector [0 ] > guess1_vector [0 ] ? guess_vector [0 ] - guess1_vector [0 ] : guess1_vector [0 ] - guess_vector [0 ];
257
+ size_t diff = guess_vector [guess_offset ] > guess1_vector [guess_offset ]
258
+ ? guess_vector [guess_offset ] - guess1_vector [guess_offset ]
259
+ : guess1_vector [guess_offset ] - guess_vector [guess_offset ];
241
260
if (diff <= 1 ) {
242
261
bool is_same = true;
243
- for (size_t i = 1 ; i < guess_arr_size - 1 ; i ++ ) {
262
+ for (size_t i = guess_offset + 1 ; i < guess_arr_size - 1 ; i ++ ) {
244
263
if (guess_vector [i ] != guess1_vector [i ]) {
245
264
is_same = false;
246
265
break ;
247
266
}
248
267
}
249
- done = is_same ;
268
+ if (is_same ) {
269
+ if (guess_precision < guess_arr_size - 1 ) {
270
+ /* If the precision has not yet reached the maximum number of digits, it will be increased. */
271
+ guess_precision = MIN (guess_precision * 3 , guess_arr_size - 1 );
272
+ updated_precision = true;
273
+ } else {
274
+ done = is_same ;
275
+ }
276
+ }
250
277
}
251
278
} while (!done );
252
279
0 commit comments