1 /*
2 Copyright (C) 2012 Fredrik Johansson
3
4 This file is part of Arb.
5
6 Arb is free software: you can redistribute it and/or modify it under
7 the terms of the GNU Lesser General Public License (LGPL) as published
8 by the Free Software Foundation; either version 2.1 of the License, or
9 (at your option) any later version. See <http://www.gnu.org/licenses/>.
10 */
11
12 #include "arb.h"
13
14 void
arb_sqrtpos(arb_t z,const arb_t x,slong prec)15 arb_sqrtpos(arb_t z, const arb_t x, slong prec)
16 {
17 if (!arb_is_finite(x))
18 {
19 if (mag_is_zero(arb_radref(x)) && arf_is_pos_inf(arb_midref(x)))
20 arb_pos_inf(z);
21 else
22 arb_zero_pm_inf(z);
23 }
24 else if (arb_contains_nonpositive(x))
25 {
26 arf_t t;
27
28 arf_init(t);
29
30 arf_set_mag(t, arb_radref(x));
31 arf_add(t, arb_midref(x), t, MAG_BITS, ARF_RND_CEIL);
32
33 if (arf_sgn(t) <= 0)
34 {
35 arb_zero(z);
36 }
37 else
38 {
39 arf_sqrt(t, t, MAG_BITS, ARF_RND_CEIL);
40 arf_mul_2exp_si(t, t, -1);
41 arf_set(arb_midref(z), t);
42 arf_get_mag(arb_radref(z), t);
43 }
44
45 arf_clear(t);
46 }
47 else
48 {
49 arb_sqrt(z, x, prec);
50 }
51
52 arb_nonnegative_part(z, z);
53 }
54
55