def factor_out_p(val, p):
    val = ZZ(val)
    if val == 0 or val == -1:
        return str(val)
    if val==1:
        return '+1'
    s = 1
    if val<0:
        s = -1
        val = -val
    ord = val.valuation(p)
    val = val/p**ord
    out = ''
    if s == -1:
        out += '-'
    else:
        out += '+'
    if ord==1:
        out +=  str(p)
    elif ord>1:
        out +=  '%d^{%d}' % (p, ord)
    if val>1:
        if ord ==1:
            out += r'\cdot '
        out += str(val)
    return out

def nu1_mu_info(w,n):
    """ return data for ST group N(U(1)) x mu(n) (of any wt > 0); these groups are not stored in the database """
    assert w > 0 and n > 0
    n = ZZ(n)
    rec = {}
    rec['label'] = "%d.2.1.d%d"%(w,n)
    rec['weight'] = w
    rec['degree'] = 2
    rec['rational'] = boolean_name(True if n <= 2 else False)
    rec['name'] = 'U(1)[D%d]'%n if n > 1 else 'N(U(1))'
    rec['pretty'] = r'\mathrm{U}(1)[D_{%d}]'%n if n > 1 else r'N(\mathrm{U}(1))'
    rec['real_dimension'] = 1
    rec['components'] = int(2*n)
    rec['ambient'] = '\mathrm{U}(2)'
    rec['connected'] = boolean_name(rec['components'] == 1)
    rec['st0_name'] = 'U(1)'
    rec['identity_component'] = st0_pretty(rec['st0_name'])
    rec['st0_description'] = '\\left\\{\\begin{bmatrix}\\alpha&0\\\\0&\\bar\\alpha\\end{bmatrix}:\\alpha\\bar\\alpha = 1,\\ \\alpha\\in\\mathbb{C}\\right\\}'
    rec['component_group'] = 'D_{%d}'%n
    rec['abelian'] = boolean_name(n <= 2)
    rec['cyclic'] = boolean_name(n <= 1)
    rec['solvable'] = boolean_name(True)
    rec['trace_zero_density']='1/2'
    rec['gens'] = r'\left\{\begin{bmatrix} 0 & 1\\ -1 & 0\end{bmatrix}, \begin{bmatrix} 1 & 0 \\ 0 & \zeta_{%d}\end{bmatrix}\right\}'%n
    rec['numgens'] = 2
    rec['subgroups'] = comma_separated_list([st_link("%d.2.1.d%d"%(w,n/p)) for p in n.prime_factors()])
    rec['supgroups'] = comma_separated_list([st_link("%d.2.1.d%d"%(w,p*n)) for p in [2,3,5]] + ["$\ldots$"])
    rec['moments'] = [['x'] + [ '\\mathrm{E}[x^{%d}]'%m for m in range(13)]]
    nu1moments = ['1','0','1','0','3','0','10','0','35','0','126','0','462']
    rec['moments'] += [['a_1'] + [nu1moments[m] if m % n == 0  else '0' for m in range(13)]]
    rec['trace_moments'] = trace_moments(rec['moments'])
    rec['counts'] = [['a_1', [[0,n]]]]
    return rec

def niceideals(F, ideals): #HNF + sage ideal + label
    """Convert a list of ideas from strongs to actual NumberField ideals

    F is a Sage NumberField

    ideals is a list of strings representing ideals I in the field, of
    the form [N,a,alpha] where N is the norm of I, a the least
    positive integer in I, and alpha a field element such that I is
    generated by a and alpha.

    The output is a list

    """
    nideals = []
    ilabel = 1
    norm = ZZ(0)
    for i in range(len(ideals)):
        N,n,idl,_ = str2ideal(F,ideals[i])
        assert idl.norm() == N and idl.smallest_integer() == n
        if N != norm:
            ilabel = ZZ(1)
            norm = N
        label = N.str() + '.' + ilabel.str()
        hnf = idl.pari_hnf().python()
        nideals.append([hnf, idl, label])
        ilabel += 1
    return nideals

def aplist(E, B=100):
    """
    Compute aplist for an elliptic curve E over Q(sqrt(5)), as a
    string->number dictionary.

    INPUT:
        - E -- an elliptic curve
        - B -- a positive integer (default: 100)

    OUTPUT:
        - dictionary mapping strings (labeled primes) to Python ints,
          with keys the primes of P with norm up to B such that the
          norm of the conductor is coprime to the characteristic of P.
    """
    from psage.modform.hilbert.sqrt5.tables import canonical_gen    
    v = {}
    from psage.modform.hilbert.sqrt5.sqrt5 import F
    labels, primes = labeled_primes_of_bounded_norm(F, B)

    from sage.all import ZZ
    N = E.conductor()
    try:
        N = ZZ(N.norm())
    except:
        N = ZZ(N)
    
    for i in range(len(primes)):
        p = primes[i]
        k = p.residue_field()
        if N.gcd(k.cardinality()) == 1:
            v[labels[i]] = int(k.cardinality() + 1 - E.change_ring(k).cardinality())
    return v

 def randomize(self, prime):
     assert not self.randomized
     prev = None
     for i in xrange(0, len(self.bp)):
         d_i_minus_one = self.bp[i].zero.nrows()
         d_i = self.bp[i].zero.ncols()
         MSZp = MatrixSpace(ZZ.residue_field(ZZ.ideal(prime)), d_i_minus_one, d_i)
         MSZp_square = MatrixSpace(ZZ.residue_field(ZZ.ideal(prime)), d_i, d_i)
         if i != 0:
             MSZp = MatrixSpace(ZZ.residue_field(ZZ.ideal(prime)), d_i_minus_one, d_i)
             self.bp[i] = self.bp[i].group(MSZp, prime).mult_left(prev.adjoint())
         if i != len(self.bp) - 1:
             cur = MSZp_square.random_element()
             self.bp[i] = self.bp[i].group(MSZp, prime).mult_right(cur)
             prev = cur
     # compute S * B_0
     d_0 = self.bp[0].zero.nrows()
     d_1 = self.bp[0].zero.ncols()
     S = matrix.identity(d_0)
     for i in xrange(d_0):
         S[i, i] = random.randint(0, prime - 1)
     MSZp = MatrixSpace(ZZ.residue_field(ZZ.ideal(prime)), d_0, d_1)
     self.bp[0] = self.bp[0].group(MSZp, prime).mult_left(S)
     # compute B_ell * T
     r = self.bp[-1].zero.nrows()
     c = self.bp[-1].zero.ncols()
     T = matrix.identity(c)
     for i in xrange(c):
         T[i, i] = random.randint(0, prime - 1)
     MSZp = MatrixSpace(ZZ.residue_field(ZZ.ideal(prime)), r, c)
     self.bp[-1] = self.bp[-1].group(MSZp, prime).mult_right(T)
     self.randomized = True

def su2_mu_info(w,n):
    """ return data for ST group SU(2) x mu(n) (of any wt > 0); these groups are not stored in the database """
    assert w > 0 and n > 0
    n = ZZ(n)
    rec = {}
    rec['label'] = "%d.2.3.c%d"%(w,n)
    rec['weight'] = w
    rec['degree'] = 2
    rec['rational'] = boolean_name(True if n <= 2 else False)
    rec['name'] = 'SU(2)[C%d]'%n if n > 1 else 'SU(2)'
    rec['pretty'] = r'\mathrm{SU}(2)[C_{%d}]'%n if n > 1 else r'\mathrm{SU}(2)'
    rec['real_dimension'] = 3
    rec['components'] = int(n)
    rec['ambient'] = '\mathrm{U}(2)'
    rec['connected'] = boolean_name(rec['components'] == 1)
    rec['st0_name'] = 'SU(2)'
    rec['identity_component'] = st0_pretty(rec['st0_name'])
    rec['st0_description'] = r'\left\{\begin{bmatrix}\alpha&\beta\\-\bar\beta&\bar\alpha\end{bmatrix}:\alpha\bar\alpha+\beta\bar\beta = 1,\ \alpha,\beta\in\mathbb{C}\right\}'
    rec['component_group'] = 'C_{%d}'%n
    rec['abelian'] = boolean_name(True)
    rec['cyclic'] = boolean_name(True)
    rec['solvable'] = boolean_name(True)
    rec['trace_zero_density']='0'
    rec['gens'] = r'\begin{bmatrix} 1 & 0 \\ 0 & \zeta_{%d}\end{bmatrix}'%n
    rec['numgens'] = 1
    rec['subgroups'] = comma_separated_list([st_link("%d.2.3.c%d"%(w,n/p)) for p in n.prime_factors()])
    rec['supgroups'] = comma_separated_list([st_link("%d.2.3.c%d"%(w,p*n)) for p in [2,3,5]] + ["$\ldots$"])
    rec['moments'] = [['x'] + [ '\\mathrm{E}[x^{%d}]'%m for m in range(13)]]
    su2moments = ['1','0','1','0','2','0','5','0','14','0','42','0','132']
    rec['moments'] += [['a_1'] + [su2moments[m] if m % n == 0  else '0' for m in range(13)]]
    rec['trace_moments'] = trace_moments(rec['moments'])
    rec['counts'] = []
    return rec

def nextchr(c):
    s = ord('a') - 1  # really 96
    l = [ZZ(ord(j) - s) for j in list(c)]
    l.reverse()
    tot = ZZ(l, 27) + 1
    newl = tot.digits(27)
    newl.reverse()
    newl = map(lambda x: x + 1 if x == 0 else x, newl)
    return ''.join([chr(j + s) for j in newl])

def render_hgm_webpage(args):
    data = None
    info = {}
    if 'label' in args:
        label = clean_input(args['label'])
        C = base.getDBConnection()
        data = C.hgm.motives.find_one({'label': label})
        if data is None:
            bread = get_bread([("Search error", url_for('.search'))])
            info['err'] = "Motive " + label + " was not found in the database."
            info['label'] = label
            return search_input_error(info, bread)
        title = 'Hypergeometric Motive:' + label
        A = data['A']
        B = data['B']
        tn,td = data['t']
        t = latex(QQ(str(tn)+'/'+str(td)))
        primes = [2, 3, 5, 7, 11, 13, 17, 19, 23, 29, 31, 37, 41, 43, 47, 53, 59, 61, 67, 71]
        locinfo = data['locinfo']
        for j in range(len(locinfo)):
            locinfo[j] = [primes[j]] + locinfo[j]
            locinfo[j][2] = poly_with_factored_coeffs(locinfo[j][2], primes[j])
        hodge = data['hodge']
        prop2 = [
            ('Degree', '\(%s\)' % data['degree']),
            ('Weight',  '\(%s\)' % data['weight']),
            ('Conductor', '\(%s\)' % data['cond']),
        ]
        # Now add factorization of conductor
        Cond = ZZ(data['cond'])
        if not (Cond.abs().is_prime() or Cond == 1):
            data['cond'] = "%s=%s" % (str(Cond), factorint(data['cond']))

        info.update({
                    'A': A,
                    'B': B,
                    't': t,
                    'degree': data['degree'],
                    'weight': data['weight'],
                    'sign': data['sign'],
                    'sig': data['sig'],
                    'hodge': hodge,
                    'cond': data['cond'],
                    'req': data['req'],
                    'locinfo': locinfo
                    })
        AB_data, t_data = data["label"].split("_t")
        #AB_data = data["label"].split("_t")[0]
        friends = [("Motive family "+AB_data.replace("_"," "), url_for(".by_family_label", label = AB_data))]
        friends.append(('L-function', url_for("l_functions.l_function_hgm_page", label=AB_data, t='t'+t_data)))
#        if rffriend != '':
#            friends.append(('Discriminant root field', rffriend))


        bread = get_bread([(label, ' ')])
        return render_template("hgm-show-motive.html", credit=HGM_credit, title=title, bread=bread, info=info, properties2=prop2, friends=friends)

def _i_func(q):
    '''Return i(B) in Katsurada's paper.
    '''
    m = ZZ(2) * (q.matrix()) ** (-1)
    i = valuation(gcd(m.list()), ZZ(2))
    m = ZZ(2) ** (-i) * m
    if all(m[a, a] % 2 == 0 for a in range(m.ncols())):
        return - i - 1
    else:
        return - i

def read_line(line, debug=0):
    r""" Parses one line from input file.  Returns and a dict containing
    fields with keys as above.

    Sample line: 11 a 1 0,-1,1,-10,-20 7 1,0 0,1,0 0,0 0,1

    Fields: label (3 fields)
            a-invariants
            p0
            For each bad  prime:  'a'                if additive
                                  lambda,mu          if multiplicative (or 'o?' if unknown)
            For each good prime:  lambda,mu          if ordinary (or 'o?' if unknown)
                                  lambda+,lambda-,mu if supersingular (or 's?' if unknown)
    """
    data = {}
    if debug: print("Parsing input line {}".format(line[:-1]))
    fields = line.split()
    label = fields[0]+fields[1]+fields[2]
    data['label'] = label
    N = ZZ(fields[0])
    badp = N.support()
    nbadp = len(badp)
    p0 = int(fields[4])
    data['iwp0'] = p0
    if debug: print("p0={}".format(p0))

    iwdata = {}

    # read data for bad primes

    for p,pdat in zip(badp,fields[5:5+nbadp]):
        p = str(p)
        if debug>1: print("p={}, pdat={}".format(p,pdat))
        if pdat in ['o?','a']:
            iwdata[p]=pdat
        else:
            iwdata[p]=[int(x) for x in pdat.split(",")]

    # read data for all primes

    for p,pdat in zip(primes(1000),fields[5+nbadp:]):
        p = str(p)
        if debug>1: print("p={}, pdat={}".format(p,pdat))
        if pdat in ['s?','o?','a']:
            iwdata[p]=pdat
        else:
            iwdata[p]=[int(x) for x in pdat.split(",")]

    data['iwdata'] = iwdata
    if debug: print("label {}, data {}".format(label,data))
    return label, data

 def includes_composite(s):
     s = s.replace(" ", "").replace("..", "-")
     for interval in s.split(","):
         if "-" in interval[1:]:
             ix = interval.index("-", 1)
             a, b = int(interval[:ix]), int(interval[ix + 1 :])
             if b == a:
                 if a != 1 and not a.is_prime():
                     return True
             if b > a and b > 3:
                 return True
         else:
             a = ZZ(interval)
             if a != 1 and not a.is_prime():
                 return True

def finite_field():
    """
    Create a random finite field with degree at most 20 and prime at most 10^6.

    OUTPUT:
        a finite field

    EXAMPLES:
        sage: sage.rings.tests.finite_field()
        Finite Field in a of size 161123^4
    """
    from sage.all import ZZ, GF
    p = ZZ.random_element(x=2, y=10**6-18).next_prime()
    d = ZZ.random_element(x=1, y=20)
    return GF(p**d,'a')

def rsa(bits):
    # only prove correctness up to 1024bits
    proof = (bits <= 1024)
    p = next_prime(ZZ.random_element(2**(bits // 2 + 1)),
                   proof=proof)
    q = next_prime(ZZ.random_element(2**(bits // 2 + 1)),
                   proof=proof)
    n = p * q
    phi_n = (p - 1) * (q - 1)
    while True:
        e = ZZ.random_element(1, phi_n)
        if gcd(e, phi_n) == 1:
            break
    d = lift(Mod(e, phi_n)**(-1))
    return e, d, n

 def includes_composite(s):
     s = s.replace(' ','').replace('..','-')
     for interval in s.split(','):
         if '-' in interval[1:]:
             ix = interval.index('-',1)
             a,b = int(interval[:ix]), int(interval[ix+1:])
             if b == a:
                 if a != 1 and not a.is_prime():
                     return True
             if b > a and b > 3:
                 return True
         else:
             a = ZZ(interval)
             if a != 1 and not a.is_prime():
                 return True

    def random_element(self, x=None, y=None, distribution=None):
        """
        Return a random integer in Pari.

        NOTE:

        The given arguments are passed to ``ZZ.random_element(...)``.

        INPUT:

        - `x`, `y` -- optional integers, that are lower and upper bound
          for the result. If only `x` is provided, then the result is
          between 0 and `x-1`, inclusive. If both are provided, then the
          result is between `x` and `y-1`, inclusive.

        - `distribution` -- optional string, so that ``ZZ`` can make sense
          of it as a probability distribution.

        EXAMPLE::

            sage: R = PariRing()
            sage: R.random_element()
            -8
            sage: R.random_element(5,13)
            12
            sage: [R.random_element(distribution="1/n") for _ in range(10)]
            [0, 1, -1, 2, 1, -95, -1, -2, -12, 0]

        """
        from sage.all import ZZ
        return self(ZZ.random_element(x,y,distribution))

def rings1():
    """
    Return an iterator over random rings.
    Return a list of pairs (f, desc), where f is a function that
    outputs a random ring that takes a ring and possibly
    some other data as constructor.

    RINGS::
        - polynomial ring in one variable over a rings0() ring.
        - polynomial ring over a rings1() ring.
        - multivariate polynomials

    EXAMPLES::

        sage: import sage.rings.tests
        sage: type(sage.rings.tests.rings0())
        <type 'list'>
    """
    v = rings0()
    X = random_rings(level=0)
    from sage.all import PolynomialRing, ZZ
    v = [(lambda : PolynomialRing(next(X), names='x'), 'univariate polynomial ring over level 0 ring'),
         (lambda : PolynomialRing(next(X), abs(ZZ.random_element(x=2,y=10)), names='x'),
                     'multivariate polynomial ring in between 2 and 10 variables over a level 0 ring')]
    return v

def relative_number_field(n=2, maxdeg=2):
    """
    Return a tower of at most n extensions each of degree at most maxdeg.

    EXAMPLES:
        sage: sage.rings.tests.relative_number_field(3)
        Number Field in aaa with defining polynomial x^2 - 15*x + 17 over its base field
    """
    from sage.all import ZZ
    K = absolute_number_field(maxdeg)
    n -= 1
    var = 'aa'
    R = ZZ['x']
    while n >= 1:
        while True:
            f = R.random_element(degree=ZZ.random_element(x=1,y=maxdeg),x=-100,y=100)
            if f.degree() <= 0: continue
            f = f * f.denominator()  # bug trac #4781
            f = f + R.gen()**maxdeg  # make monic
            if f.is_irreducible():
                break
        K = K.extension(f,var)
        var += 'a'
        n -= 1
    return K

def allbsd(line):
    r""" Parses one line from an allbsd file.  Returns the label and a
    dict containing fields with keys 'conductor', 'iso', 'number',
    'ainvs', 'rank', 'torsion', 'torsion_primes', 'tamagawa_product',
    'real_period', 'special_value', 'regulator', 'sha_an', 'sha',
    'sha_primes', all values being strings or floats or ints or lists
    of ints.

    Input line fields:

    conductor iso number ainvs rank torsion tamagawa_product real_period special_value regulator sha_an

    Sample input line:

    11 a 1 [0,-1,1,-10,-20] 0 5 5 1.2692093042795534217 0.25384186085591068434 1 1.00000000000000000000

    """
    data = split(line)
    label = data[0] + data[1] + data[2]
    ainvs = parse_ainvs(data[3])

    torsion = ZZ(data[5])
    sha_an = RR(data[10])
    sha = sha_an.round()
    sha_primes = sha.prime_divisors()
    torsion_primes = torsion.prime_divisors()

    data = {
        'conductor': int(data[0]),
        'iso': data[0] + data[1],
        'number': int(data[2]),
        'ainvs': ainvs,
        'rank': int(data[4]),
        'tamagawa_product': int(data[6]),
        'real_period': float(data[7]),
        'special_value': float(data[8]),
        'regulator': float(data[9]),
        'sha_an': float(sha_an),
        'sha':  int(sha),
        'sha_primes':  [int(p) for p in sha_primes],
        'torsion':  int(torsion),
        'torsion_primes':  [int(p) for p in torsion_primes]
        }

    return label, data

    def init_dir_char(self, chi):
        """
        Initiate with a Web Dirichlet character.
        """
        self.original_object = [chi]
        chi = chi.chi.primitive_character()
        self.object_type = "dirichletcharacter"
        self.dim = 1
        self.motivic_weight = 0
        self.conductor = ZZ(chi.conductor())
        self.bad_semistable_primes = []
        self.bad_pot_good = self.conductor.prime_factors()
        if chi.is_odd():
            aa = 1
            bb = I
        else:
            aa = 0
            bb = 1
        self.sign = chi.gauss_sum_numerical() / (bb * float(sqrt(chi.modulus())) )
        # this has now type python complex. later we need a gp complex
        self.sign = ComplexField()(self.sign)
        self.mu_fe = [aa]
        self.nu_fe = []
        self.gammaV = [aa]
        self.langlands = True
        self.selfdual = (chi.multiplicative_order() <= 2)
        # rather than all(  abs(chi(m).imag) < 0.0001 for m in range(chi.modulus() ) )
        self.primitive = True
        self.set_dokchitser_Lfunction()
        self.set_number_of_coefficients()
        self.dirichlet_coefficients = [ chi(m) for m in range(self.numcoeff + 1) ]
        if self.selfdual:
            self.coefficient_type = 2
        else:
            self.coefficient_type = 3
        self.coefficient_period = chi.modulus()
        self.besancon_bound = 10000

        def eu(p):
            """
            local euler factor
            """
            if self.selfdual:
                K = QQ
            else:
                K = ComplexField()
            R = PolynomialRing(K, "T")
            T = R.gens()[0]
            if self.conductor % p != 0:
                return  1 - ComplexField()(chi(p)) * T
            else:
                return R(1)

        self.local_euler_factor = eu
        self.ld.gp().quit()

def twoadic(line):
    r""" Parses one line from a 2adic file.  Returns the label and a dict
    containing fields with keys '2adic_index', '2adic_log_level',
    '2adic_gens' and '2adic_label'.

    Input line fields:

    conductor iso number ainvs index level gens label

    Sample input lines:

    110005 a 2 [1,-1,1,-185793,29503856] 12 4 [[3,0,0,1],[3,2,2,3],[3,0,0,3]] X24
    27 a 1 [0,0,1,0,-7] inf inf [] CM
    """
    data = split(line)
    assert len(data)==8
    label = data[0] + data[1] + data[2]
    model = data[7]
    if model == 'CM':
        return label, {
            '2adic_index': int(0),
            '2adic_log_level': None,
            '2adic_gens': None,
            '2adic_label': None,
        }

    index = int(data[4])
    level = ZZ(data[5])
    log_level = int(level.valuation(2))
    assert 2**log_level==level
    if data[6]=='[]':
        gens=[]
    else:
        gens = data[6][1:-1].replace('],[','];[').split(';')
        gens = [[int(c) for c in g[1:-1].split(',')] for g in gens]

    return label, {
            '2adic_index': index,
            '2adic_log_level': log_level,
            '2adic_gens': gens,
            '2adic_label': model,
    }