# Ticket #10720(needs_work enhancement)

Opened 2 years ago

## nth_root in power series

Reported by: Owned by: pernici pernici minor sage-5.10 commutative algebra power series robertwb N/A mario pernici

computation of an nth root of a power series using the Newton method

Apply

## Change History

### comment:1 Changed 2 years ago by pernici

• Description modified (diff)
• Authors set to mario pernici
• Component changed from PLEASE CHANGE to commutative algebra
• Priority changed from major to minor
• Owner changed from tbd to pernici
• Milestone set to sage-4.6.2

### comment:2 Changed 2 years ago by pernici

• Description modified (diff)

### comment:3 Changed 2 years ago by pernici

The nth-root of power series y = x^n is computed using the Newton method for y = x^-n

```x' = (1+1/n)*x - y*x^(n+1)/n
```
```sage: R.<t> = QQ[]
sage: p = (1 + 2*t + 5*t^2 + 7*t^3 + O(t^4))^3
sage: p.nth_root(3)
1 + 2*t + 5*t^2 + 7*t^3 + O(t^4)
sage: p = (1 + 2*t + 5*t^2 + 7*t^3 + O(t^4))^-3
sage: p.nth_root(-3)
1 + 2*t + 5*t^2 + 7*t^3 + O(t^4)
```

With this patch nth_root works also for series on ZZ

```sage: R.<t> = ZZ[[]]
sage: p = 1 + 4*t^2 + 3*t^3 + O(t^4)
sage: p.nth_root(2)
1 + 2*t^2 + 3/2*t^3 + O(t^4)
sage: p.sqrt(2)
1 + 2*t^2 + 3/2*t^3 + O(t^4)
sage: p.nth_root(3)
1 + 4/3*t^2 + t^3 + O(t^4)
sage: p = 4*t^2 + 3*t^3 + O(t^4)
sage: p.nth_root(2)
2*t + 3/4*t^2 + O(t^3)
sage: p.sqrt(2)
Traceback (most recent call last):
...
TypeError: no conversion of this rational to integer
```

there is a bug in sqrt in this case.

nth_root can be used to compute the square root of series which currently sqrt does not support

```sage: R.<x,y> = QQ[]
sage: S.<t> = R[[]]
sage: p = 4 + t*x + t^2*y + O(t^3)
sage: p.nth_root(2)
2 + 1/4*x*t + (-1/64*x^2 + 1/4*y)*t^2 + O(t^3)
sage: p = 32*t^5 + x*t^6 + y*t^7 + O(t^8)
sage: p.nth_root(5)
2*t + 1/80*x*t^2 + (-1/6400*x^2 + 1/80*y)*t^3 + O(t^4)
sage: p.sqrt()
Traceback (most recent call last):
...
AttributeError: 'sage.rings.polynomial.multi_polynomial_libsingular.MPolynomial_libsingular' object has no attribute 'sqrt'
```

nth-root iteration formula is division-free, which is convenient in this case; in this case sqrt could call nth_root.

nth_root can be used in the multivariate series considered in ticket #1956 .

### comment:4 Changed 2 years ago by pernici

in the above message there is a wrong example on where sqrt fails; here is the corrected example

```sage: R.<x,y> = QQ[]
sage: S.<t> = R[[]]
sage: p = 4 + t*x + t^2*y + O(t^3)
sage: p.nth_root(2)
2 + 1/4*x*t + (-1/64*x^2 + 1/4*y)*t^2 + O(t^3)
sage: p.sqrt()
Traceback (most recent call last):
...
AttributeError: 'sage.rings.polynomial.multi_polynomial_libsingular.MPolynomial_libsingular' object has no attribute 'sqrt'
}
```

### comment:5 Changed 2 years ago by robertwb

• Status changed from new to needs_review
• Type changed from PLEASE CHANGE to defect

### comment:6 Changed 2 years ago by robertwb

Apply only trac_10720_power_series_nth_root_2.patch

### comment:7 Changed 2 years ago by pernici

• Description modified (diff)

With trac_10720_power_series_nth_root_3.patch nth_root is much faster for n large

```sage: S.<t> = QQ[[]]
sage: p = t.exp(30)
sage: %time p1 = p.nth_root(1000)
Wall time: 0.01 s
```

with the previous version it takes 84s

### comment:8 Changed 2 years ago by pernici

With trac_10720_power_series_nth_root_4.patch nth_root and inversion become faster for series on polynomial rings, when used together with trac_10480_fast_PowerSeries_poly_multiplication-alternative.patch from ticket #10480 (otherwise these functions perform as without the patch).

In the following benchmarks apply trac_karatsuba_improvements.patch and trac_10480_fast_PowerSeries_poly_multiplication-alternative.patch, trac_10720_power_series_nth_root_2.patch and trac_10720_power_series_nth_root_3.patch; in (2) apply also the new patch

```sage: N = 2
sage: R= PolynomialRing(QQ,N,'x')
sage: x = R.gens()
sage: S.<t> = R[[]]
sage: sx = sum(x)
sage: prec = 15
sage: p = (t.exp(prec)).subs({t:t*sx})
sage: %time p1 = p^-1
sage: %time p1 = p.nth_root(2)
```

Inversion benchmark:

```N  prec  (1)    (2)
2  15    0.06   0.01
3  15    2.1    0.07
4  15    58     0.5
5  15    1440   3.4
```

Square root benchmark:

```N  prec  (1)    (2)
2  15    0.15   0.03
3  15    5.2    0.31
4  15    161    3.3
```

The speedup is large because the terms with high degree in t are large polynomials, due to the combinatorial possibilities with several variables; using _mul_ instead of _mul_trunc several products involving these large polynomials are computed and later truncated away.

### comment:9 Changed 2 years ago by pernici

• Description modified (diff)

### comment:10 Changed 19 months ago by johanbosman

• Type changed from defect to enhancement

This is an enhancement rather than a defect.

### comment:11 Changed 17 months ago by vbraun

Which patches need to be applied? Please put some information in the ticket description.

If trac_10720_power_series_nth_root_3.patch is meant to be applied, please add a doctest to __pow_trunc.

### comment:13 Changed 3 months ago by saraedum

• Status changed from needs_review to needs_work
• Description modified (diff)

As Volker mentioned, __pow_trunc lacks a docstring. (Just to remove this from the list of tickets needing review until that is fixed.)

### comment:14 Changed 2 months ago by chapoton

One related remark:

Thanks to ticket #14115, it is now easy to implement a general rational power by

```def rational_power(f,alpha):
return (f.log()*alpha).exp()
```

I do not know how it compares in terms of speed with the method of this ticket.

Note: See TracTickets for help on using tickets.