Index: sageRunSLZ.sage
===================================================================
--- sageRunSLZ.sage	(revision 212)
+++ sageRunSLZ.sage	(revision 213)
@@ -1,4 +1,509 @@
+r"""
+Main SLZ algorithm body in several versions.
+
+AUTHORS:
+- S.T. (2015-10-10): initial version
+
+Examples:
+    TODO
 """
+print "sageRationalOperations loading..."
+
+def srs_compute_lattice_volume(inputFunction,
+                    inputLowerBound,
+                    inputUpperBound,
+                    alpha,
+                    degree,
+                    precision,
+                    emin,
+                    emax,
+                    targetHardnessToRound,
+                    debug = False):
+    """
+    Changes from V2:
+        Root search is changed:
+            - we compute the resultants in i and in t;
+            - we compute the roots set of each of these resultants;
+            - we combine all the possible pairs between the two sets;
+            - we check these pairs in polynomials for correctness.
+    Changes from V1: 
+        1- check for roots as soon as a resultant is computed;
+        2- once a non null resultant is found, check for roots;
+        3- constant resultant == no root.
+    """
+
+    if debug:
+        print "Function                :", inputFunction
+        print "Lower bound             :", inputLowerBound
+        print "Upper bounds            :", inputUpperBound
+        print "Alpha                   :", alpha
+        print "Degree                  :", degree
+        print "Precision               :", precision
+        print "Emin                    :", emin
+        print "Emax                    :", emax
+        print "Target hardness-to-round:", targetHardnessToRound
+        print
+    ## Important constants.
+    ### Stretch the interval if no error happens.
+    noErrorIntervalStretch = 1 + 2^(-5)
+    ### If no vector validates the Coppersmith condition, shrink the interval
+    #   by the following factor.
+    noCoppersmithIntervalShrink = 1/2
+    ### If only (or at least) one vector validates the Coppersmith condition, 
+    #   shrink the interval by the following factor.
+    oneCoppersmithIntervalShrink = 3/4
+    #### If only null resultants are found, shrink the interval by the 
+    #    following factor.
+    onlyNullResultantsShrink     = 3/4
+    ## Structures.
+    RRR         = RealField(precision)
+    RRIF        = RealIntervalField(precision)
+    ## Converting input bound into the "right" field.
+    lowerBound = RRR(inputLowerBound)
+    upperBound = RRR(inputUpperBound)
+    ## Before going any further, check domain and image binade conditions.
+    print inputFunction(1).n()
+    output = slz_fix_bounds_for_binades(lowerBound, upperBound, inputFunction)
+    if output is None:
+        print "Invalid domain/image binades. Domain:",\
+        lowerBound, upperBound, "Images:", \
+        inputFunction(lowerBound), inputFunction(upperBound)
+        raise Exception("Invalid domain/image binades.")
+    lb = output[0] ; ub = output[1]
+    if lb != lowerBound or ub != upperBound:
+        print "lb:", lb, " - ub:", ub
+        print "Invalid domain/image binades. Domain:",\
+        lowerBound, upperBound, "Images:", \
+        inputFunction(lowerBound), inputFunction(upperBound)
+        raise Exception("Invalid domain/image binades.")
+    #
+    ## Progam initialization
+    ### Approximation polynomial accuracy and hardness to round.
+    polyApproxAccur           = 2^(-(targetHardnessToRound + 1))
+    polyTargetHardnessToRound = targetHardnessToRound + 1
+    ### Significand to integer conversion ratio.
+    toIntegerFactor           = 2^(precision-1)
+    print "Polynomial approximation required accuracy:", polyApproxAccur.n()
+    ### Variables and rings for polynomials and root searching.
+    i=var('i')
+    t=var('t')
+    inputFunctionVariable = inputFunction.variables()[0]
+    function = inputFunction.subs({inputFunctionVariable:i})
+    #  Polynomial Rings over the integers, for root finding.
+    Zi = ZZ[i]
+    Zt = ZZ[t]
+    Zit = ZZ[i,t]
+    ## Number of iterations limit.
+    maxIter = 100000
+    #
+    ## Compute the scaled function and the degree, in their Sollya version 
+    #  once for all.
+    (scaledf, sdlb, sdub, silb, siub) = \
+        slz_compute_scaled_function(function, lowerBound, upperBound, precision)
+    print "Scaled function:", scaledf._assume_str().replace('_SAGE_VAR_', '')
+    scaledfSo = sollya_lib_parse_string(scaledf._assume_str().replace('_SAGE_VAR_', ''))
+    degreeSo  = pobyso_constant_from_int_sa_so(degree)
+    #
+    ## Compute the scaling. boundsIntervalRifSa defined out of the loops.
+    domainBoundsInterval   = RRIF(lowerBound, upperBound)
+    (unscalingFunction, scalingFunction) = \
+        slz_interval_scaling_expression(domainBoundsInterval, i)
+    #print scalingFunction, unscalingFunction
+    ## Set the Sollya internal precision (with an arbitrary minimum of 192).
+    internalSollyaPrec = ceil((RR('1.5') * targetHardnessToRound) / 64) * 64
+    if internalSollyaPrec < 192:
+        internalSollyaPrec = 192
+    pobyso_set_prec_sa_so(internalSollyaPrec)
+    print "Sollya internal precision:", internalSollyaPrec
+    ## Some variables.
+    ### General variables
+    lb                  = sdlb
+    ub                  = sdub
+    nbw                 = 0
+    intervalUlp         = ub.ulp()
+    #### Will be set by slz_interval_and_polynomila_to_sage.
+    ic                  = 0 
+    icAsInt             = 0    # Set from ic.
+    solutionsSet        = set()
+    tsErrorWidth        = []
+    csErrorVectors      = []
+    csVectorsResultants = []
+    floatP              = 0  # Taylor polynomial.
+    floatPcv            = 0  # Ditto with variable change.
+    intvl               = "" # Taylor interval
+    terr                = 0  # Taylor error.
+    iterCount = 0
+    htrnSet = set()
+    ### Timers and counters.
+    wallTimeStart                   = 0
+    cpuTimeStart                    = 0
+    taylCondFailedCount             = 0
+    coppCondFailedCount             = 0
+    resultCondFailedCount           = 0
+    coppCondFailed                  = False
+    resultCondFailed                = False
+    globalResultsList               = []
+    basisConstructionsCount         = 0
+    basisConstructionsFullTime      = 0
+    basisConstructionTime           = 0
+    reductionsCount                 = 0
+    reductionsFullTime              = 0
+    reductionTime                   = 0
+    resultantsComputationsCount     = 0
+    resultantsComputationsFullTime  = 0
+    resultantsComputationTime       = 0
+    rootsComputationsCount          = 0
+    rootsComputationsFullTime       = 0
+    rootsComputationTime            = 0
+    print
+    ## Global times are started here.
+    wallTimeStart                   = walltime()
+    cpuTimeStart                    = cputime()
+    ## Main loop.
+    while True:
+        if lb >= sdub:
+            print "Lower bound reached upper bound."
+            break
+        if iterCount == maxIter:
+            print "Reached maxIter. Aborting"
+            break
+        iterCount += 1
+        print "[", lb, ",", ub, "]", ((ub - lb) / intervalUlp).log2().n(), \
+            "log2(numbers)." 
+        ### Compute a Sollya polynomial that will honor the Taylor condition.
+        prceSo = slz_compute_polynomial_and_interval(scaledfSo,
+                                                     degreeSo,
+                                                     lb,
+                                                     ub,
+                                                     polyApproxAccur)
+        ### Convert back the data into Sage space.                         
+        (floatP, floatPcv, intvl, ic, terr) = \
+            slz_interval_and_polynomial_to_sage((prceSo[0], prceSo[0],
+                                                 prceSo[1], prceSo[2], 
+                                                 prceSo[3]))
+        intvl = RRIF(intvl)
+        ## Clean-up Sollya stuff.
+        for elem in prceSo:
+            sollya_lib_clear_obj(elem)
+        #print  floatP, floatPcv, intvl, ic, terr
+        #print floatP
+        #print intvl.endpoints()[0].n(), \
+        #      ic.n(),
+        #intvl.endpoints()[1].n()
+        ### Check returned data.
+        #### Is approximation error OK?
+        if terr > polyApproxAccur:
+            exceptionErrorMess  = \
+                "Approximation failed - computed error:" + \
+                str(terr) + " - target error: "
+            exceptionErrorMess += \
+                str(polyApproxAccur) + ". Aborting!"
+            raise Exception(exceptionErrorMess)
+        #### Is lower bound OK?
+        if lb != intvl.endpoints()[0]:
+            exceptionErrorMess =  "Wrong lower bound:" + \
+                               str(lb) + ". Aborting!"
+            raise Exception(exceptionErrorMess)
+        #### Set upper bound.
+        if ub > intvl.endpoints()[1]:
+            ub = intvl.endpoints()[1]
+            print "[", lb, ",", ub, "]", ((ub - lb) / intervalUlp).log2().n(), \
+            "log2(numbers)." 
+            taylCondFailedCount += 1
+        #### Is interval not degenerate?
+        if lb >= ub:
+            exceptionErrorMess = "Degenerate interval: " + \
+                                "lowerBound(" + str(lb) +\
+                                ")>= upperBound(" + str(ub) + \
+                                "). Aborting!"
+            raise Exception(exceptionErrorMess)
+        #### Is interval center ok?
+        if ic <= lb or ic >= ub:
+            exceptionErrorMess =  "Invalid interval center for " + \
+                                str(lb) + ',' + str(ic) + ',' +  \
+                                str(ub) + ". Aborting!"
+            raise Exception(exceptionErrorMess)
+        ##### Current interval width and reset future interval width.
+        bw = ub - lb
+        nbw = 0
+        icAsInt = int(ic * toIntegerFactor)
+        #### The following ratio is always >= 1. In case we may want to
+        #    enlarge the interval
+        curTaylErrRat = polyApproxAccur / terr
+        ### Make the  integral transformations.
+        #### Bounds and interval center.
+        intIc = int(ic * toIntegerFactor)
+        intLb = int(lb * toIntegerFactor) - intIc
+        intUb = int(ub * toIntegerFactor) - intIc
+        #
+        #### Polynomials
+        basisConstructionTime         = cputime()
+        ##### To a polynomial with rational coefficients with rational arguments
+        ratRatP = slz_float_poly_of_float_to_rat_poly_of_rat_pow_two(floatP)
+        ##### To a polynomial with rational coefficients with integer arguments
+        ratIntP = \
+            slz_rat_poly_of_rat_to_rat_poly_of_int(ratRatP, precision)
+        #####  Ultimately a multivariate polynomial with integer coefficients  
+        #      with integer arguments.
+        coppersmithTuple = \
+            slz_rat_poly_of_int_to_poly_for_coppersmith(ratIntP, 
+                                                        precision, 
+                                                        targetHardnessToRound, 
+                                                        i, t)
+        #### Recover Coppersmith information.
+        intIntP = coppersmithTuple[0]
+        N = coppersmithTuple[1]
+        nAtAlpha = N^alpha
+        tBound = coppersmithTuple[2]
+        leastCommonMultiple = coppersmithTuple[3]
+        iBound = max(abs(intLb),abs(intUb))
+        basisConstructionsFullTime        += cputime(basisConstructionTime)
+        basisConstructionsCount           += 1
+        reductionTime                     = cputime()
+        #### Compute the reduced polynomials.
+        ccReducedPolynomialsList =  \
+            slz_compute_coppersmith_reduced_polynomials_with_lattice_volume(intIntP, 
+                                                        alpha, 
+                                                        N, 
+                                                        iBound, 
+                                                        tBound)
+        if ccReducedPolynomialsList is None:
+            raise Exception("Reduction failed.")
+        reductionsFullTime    += cputime(reductionTime)
+        reductionsCount       += 1
+        if len(ccReducedPolynomialsList) < 2:
+            print "Nothing to form resultants with."
+            print
+            coppCondFailedCount += 1
+            coppCondFailed = True
+            ##### Apply a different shrink factor according to 
+            #  the number of compliant polynomials.
+            if len(ccReducedPolynomialsList) == 0:
+                ub = lb + bw * noCoppersmithIntervalShrink
+            else: # At least one compliant polynomial.
+                ub = lb + bw * oneCoppersmithIntervalShrink
+            if ub > sdub:
+                ub = sdub
+            if lb == ub:
+                raise Exception("Cant shrink interval \
+                anymore to get Coppersmith condition.")
+            nbw = 0
+            continue
+        #### We have at least two polynomials.
+        #    Let us try to compute resultants.
+        #    For each resultant computed, go for the solutions.
+        ##### Build the pairs list.
+        polyPairsList           = []
+        for polyOuterIndex in xrange(0, len(ccReducedPolynomialsList) - 1):
+            for polyInnerIndex in xrange(polyOuterIndex+1, 
+                                         len(ccReducedPolynomialsList)):
+                polyPairsList.append((ccReducedPolynomialsList[polyOuterIndex],
+                                      ccReducedPolynomialsList[polyInnerIndex]))
+        #### Actual root search.
+        rootsSet            = set()
+        hasNonNullResultant = False
+        for polyPair in polyPairsList:
+            if hasNonNullResultant:
+                break
+            resultantsComputationTime   = cputime()
+            currentResultantI = \
+                slz_resultant(polyPair[0],
+                              polyPair[1],
+                              t)
+            resultantsComputationsCount    += 1
+            if currentResultantI is None:
+                resultantsComputationsFullTime +=  \
+                    cputime(resultantsComputationTime)
+                print "Nul resultant"
+                continue # Next polyPair.
+            currentResultantT = \
+                slz_resultant(polyPair[0],
+                              polyPair[1],
+                              i)
+            resultantsComputationsFullTime += cputime(resultantsComputationTime)
+            resultantsComputationsCount    += 1
+            if currentResultantT is None:                    
+                print "Nul resultant"
+                continue # Next polyPair.
+            else:
+                hasNonNullResultant = True
+            #### We have a non null resultants pair. From now on, whatever the
+            #    root search yields, no extra root search is necessary.
+            #### A constant resultant leads to no root. Root search is done.
+            if currentResultantI.degree() < 1:
+                print "Resultant is constant:", currentResultantI
+                break # Next polyPair and should break.
+            if currentResultantT.degree() < 1:
+                print "Resultant is constant:", currentResultantT
+                break # Next polyPair and should break.
+            #### Actual roots computation.
+            rootsComputationTime       = cputime()
+            ##### Compute i roots
+            iRootsList = Zi(currentResultantI).roots()
+            rootsComputationsCount      +=  1
+            if len(iRootsList) == 0:
+                rootsComputationsFullTime   =   cputime(rootsComputationTime)
+                print "No roots in \"i\"."
+                break # No roots in i.
+            tRootsList = Zt(currentResultantT).roots()
+            rootsComputationsFullTime   =   cputime(rootsComputationTime)
+            rootsComputationsCount      +=  1
+            if len(tRootsList) == 0:
+                print "No roots in \"t\"."
+                break # No roots in i.
+            ##### For each iRoot, get a tRoot and check against the polynomials.
+            for iRoot in iRootsList:
+                ####### Roots returned by roots() are (value, multiplicity) 
+                #       tuples.
+                #print "iRoot:", iRoot
+                for tRoot in tRootsList:
+                ###### Use the tRoot against each polynomial, alternatively.
+                    if polyPair[0](iRoot[0],tRoot[0]) != 0:
+                        continue
+                    if polyPair[1](iRoot[0],tRoot[0]) != 0:
+                        continue
+                    rootsSet.add((iRoot[0], tRoot[0]))
+            # End of roots computation.
+        # End loop for polyPair in polyParsList.  Will break at next iteration.
+        # since a non null resultant was found.
+        #### Prepare for results for the current interval..
+        intervalResultsList = []
+        intervalResultsList.append((lb, ub))
+        #### Check roots.
+        rootsResultsList = []
+        for root in rootsSet:
+            specificRootResultsList = []
+            failingBounds = []
+            intIntPdivN = intIntP(root[0], root[1]) / N
+            if int(intIntPdivN) != intIntPdivN:
+                continue # Next root
+            # Root qualifies for modular equation, test it for hardness to round.
+            hardToRoundCaseAsFloat = RRR((icAsInt + root[0]) / toIntegerFactor)
+            #print "Before unscaling:", hardToRoundCaseAsFloat.n(prec=precision)
+            #print scalingFunction
+            scaledHardToRoundCaseAsFloat = \
+                scalingFunction(hardToRoundCaseAsFloat) 
+            print "Candidate HTRNc at x =", \
+                scaledHardToRoundCaseAsFloat.n().str(base=2),
+            if slz_is_htrn(scaledHardToRoundCaseAsFloat,
+                           function,
+                           2^-(targetHardnessToRound),
+                           RRR):
+                print hardToRoundCaseAsFloat, "is HTRN case."
+                if lb <= hardToRoundCaseAsFloat and hardToRoundCaseAsFloat <= ub:
+                    print "Found in interval."
+                else:
+                    print "Found out of interval."
+                specificRootResultsList.append(hardToRoundCaseAsFloat.n().str(base=2))
+                # Check the root is in the bounds
+                if abs(root[0]) > iBound or abs(root[1]) > tBound:
+                    print "Root", root, "is out of bounds for modular equation."
+                    if abs(root[0]) > iBound:
+                        print "root[0]:", root[0]
+                        print "i bound:", iBound
+                        failingBounds.append('i')
+                        failingBounds.append(root[0])
+                        failingBounds.append(iBound)
+                    if abs(root[1]) > tBound:
+                        print "root[1]:", root[1]
+                        print "t bound:", tBound
+                        failingBounds.append('t')
+                        failingBounds.append(root[1])
+                        failingBounds.append(tBound)
+                if len(failingBounds) > 0:
+                    specificRootResultsList.append(failingBounds)
+            else: # From slz_is_htrn...
+                print "is not an HTRN case."
+            if len(specificRootResultsList) > 0:
+                rootsResultsList.append(specificRootResultsList)
+        if len(rootsResultsList) > 0:
+            intervalResultsList.append(rootsResultsList)
+        ### Check if a non null resultant was found. If not shrink the interval.
+        if not hasNonNullResultant:
+            print "Only null resultants for this reduction, shrinking interval."
+            resultCondFailed      =  True
+            resultCondFailedCount += 1
+            ### Shrink interval for next iteration.
+            ub = lb + bw * onlyNullResultantsShrink
+            if ub > sdub:
+                ub = sdub
+            nbw = 0
+            continue
+        #### An intervalResultsList has at least the bounds.
+        globalResultsList.append(intervalResultsList)   
+        #### Compute an incremented width for next upper bound, only
+        #    if not Coppersmith condition nor resultant condition
+        #    failed at the previous run. 
+        if not coppCondFailed and not resultCondFailed:
+            nbw = noErrorIntervalStretch * bw
+        else:
+            nbw = bw
+        ##### Reset the failure flags. They will be raised
+        #     again if needed.
+        coppCondFailed   = False
+        resultCondFailed = False
+        #### For next iteration (at end of loop)
+        #print "nbw:", nbw
+        lb  = ub
+        ub += nbw     
+        if ub > sdub:
+            ub = sdub
+        print
+    # End while True
+    ## Main loop just ended.
+    globalWallTime = walltime(wallTimeStart)
+    globalCpuTime  = cputime(cpuTimeStart)
+    ## Output results
+    print ; print "Intervals and HTRNs" ; print
+    for intervalResultsList in globalResultsList:
+        print "[", intervalResultsList[0][0], ",",intervalResultsList[0][1], "]",
+        if len(intervalResultsList) > 1:
+            rootsResultsList = intervalResultsList[1]
+            for specificRootResultsList in rootsResultsList:
+                print "\t", specificRootResultsList[0],
+                if len(specificRootResultsList) > 1:
+                    print specificRootResultsList[1],
+        print ; print
+    #print globalResultsList
+    #
+    print "Timers and counters"
+    print
+    print "Number of iterations:", iterCount
+    print "Taylor condition failures:", taylCondFailedCount
+    print "Coppersmith condition failures:", coppCondFailedCount
+    print "Resultant condition failures:", resultCondFailedCount
+    print "Iterations count: ", iterCount
+    print "Number of intervals:", len(globalResultsList)
+    print "Number of basis constructions:", basisConstructionsCount 
+    print "Total CPU time spent in basis constructions:", \
+        basisConstructionsFullTime
+    if basisConstructionsCount != 0:
+        print "Average basis construction CPU time:", \
+            basisConstructionsFullTime/basisConstructionsCount
+    print "Number of reductions:", reductionsCount
+    print "Total CPU time spent in reductions:", reductionsFullTime
+    if reductionsCount != 0:
+        print "Average reduction CPU time:", \
+            reductionsFullTime/reductionsCount
+    print "Number of resultants computation rounds:", \
+        resultantsComputationsCount
+    print "Total CPU time spent in resultants computation rounds:", \
+        resultantsComputationsFullTime
+    if resultantsComputationsCount != 0:
+        print "Average resultants computation round CPU time:", \
+            resultantsComputationsFullTime/resultantsComputationsCount
+    print "Number of root finding rounds:", rootsComputationsCount
+    print "Total CPU time spent in roots finding rounds:", \
+        rootsComputationsFullTime
+    if rootsComputationsCount != 0:
+        print "Average roots finding round CPU time:", \
+            rootsComputationsFullTime/rootsComputationsCount
+    print "Global Wall time:", globalWallTime
+    print "Global CPU time:", globalCpuTime
+    ## Output counters
+# End srs_compute_lattice_volume
+ 
+"""
 SLZ runtime function.
 """
 
@@ -1431,7 +1936,7 @@
     ## Output counters
 # End srs_runSLZ-v03
 
-def srs_compute_lattice_volume(inputFunction,
+def srs_run_SLZ_v04(inputFunction,
                     inputLowerBound,
                     inputUpperBound,
                     alpha,
@@ -1442,6 +1947,11 @@
                     targetHardnessToRound,
                     debug = False):
     """
+    Changes from V3:
+        Root search is changed again:
+            - only resultants in i are computed;
+            - root are searched for;
+            - if any, they are tested for hardness-to-round.
     Changes from V2:
         Root search is changed:
             - we compute the resultants in i and in t;
@@ -1684,11 +2194,11 @@
         reductionTime                     = cputime()
         #### Compute the reduced polynomials.
         ccReducedPolynomialsList =  \
-            slz_compute_coppersmith_reduced_polynomials_with_lattice_volume(intIntP, 
-                                                        alpha, 
-                                                        N, 
-                                                        iBound, 
-                                                        tBound)
+                slz_compute_coppersmith_reduced_polynomials(intIntP, 
+                                                            alpha, 
+                                                            N, 
+                                                            iBound, 
+                                                            tBound)
         if ccReducedPolynomialsList is None:
             raise Exception("Reduction failed.")
         reductionsFullTime    += cputime(reductionTime)
@@ -1722,85 +2232,54 @@
                 polyPairsList.append((ccReducedPolynomialsList[polyOuterIndex],
                                       ccReducedPolynomialsList[polyInnerIndex]))
         #### Actual root search.
-        rootsSet            = set()
+        iRootsSet           = set()
         hasNonNullResultant = False
         for polyPair in polyPairsList:
-            if hasNonNullResultant:
-                break
             resultantsComputationTime   = cputime()
             currentResultantI = \
                 slz_resultant(polyPair[0],
                               polyPair[1],
                               t)
             resultantsComputationsCount    += 1
+            resultantsComputationsFullTime +=  \
+                cputime(resultantsComputationTime)
+            #### Function slz_resultant returns None both for None and O
+            #    resultants.
             if currentResultantI is None:
-                resultantsComputationsFullTime +=  \
-                    cputime(resultantsComputationTime)
                 print "Nul resultant"
                 continue # Next polyPair.
-            currentResultantT = \
-                slz_resultant(polyPair[0],
-                              polyPair[1],
-                              i)
-            resultantsComputationsFullTime += cputime(resultantsComputationTime)
-            resultantsComputationsCount    += 1
-            if currentResultantT is None:                    
-                print "Nul resultant"
-                continue # Next polyPair.
-            else:
-                hasNonNullResultant = True
-            #### We have a non null resultants pair. From now on, whatever the
+            ## We deleted the currentResultantI computation.
+            #### We have a non null resultant. From now on, whatever this
             #    root search yields, no extra root search is necessary.
+            hasNonNullResultant = True
             #### A constant resultant leads to no root. Root search is done.
             if currentResultantI.degree() < 1:
                 print "Resultant is constant:", currentResultantI
-                break # Next polyPair and should break.
-            if currentResultantT.degree() < 1:
-                print "Resultant is constant:", currentResultantT
-                break # Next polyPair and should break.
-            #### Actual roots computation.
-            rootsComputationTime       = cputime()
-            ##### Compute i roots
+                break # There is no root.
+            #### Actual iroots computation.
+            rootsComputationTime        = cputime()
             iRootsList = Zi(currentResultantI).roots()
             rootsComputationsCount      +=  1
+            rootsComputationsFullTime   =   cputime(rootsComputationTime)
             if len(iRootsList) == 0:
-                rootsComputationsFullTime   =   cputime(rootsComputationTime)
                 print "No roots in \"i\"."
                 break # No roots in i.
-            tRootsList = Zt(currentResultantT).roots()
-            rootsComputationsFullTime   =   cputime(rootsComputationTime)
-            rootsComputationsCount      +=  1
-            if len(tRootsList) == 0:
-                print "No roots in \"t\"."
-                break # No roots in i.
-            ##### For each iRoot, get a tRoot and check against the polynomials.
-            for iRoot in iRootsList:
-                ####### Roots returned by roots() are (value, multiplicity) 
-                #       tuples.
-                #print "iRoot:", iRoot
-                for tRoot in tRootsList:
-                ###### Use the tRoot against each polynomial, alternatively.
-                    if polyPair[0](iRoot[0],tRoot[0]) != 0:
-                        continue
-                    if polyPair[1](iRoot[0],tRoot[0]) != 0:
-                        continue
-                    rootsSet.add((iRoot[0], tRoot[0]))
-            # End of roots computation.
-        # End loop for polyPair in polyParsList.  Will break at next iteration.
-        # since a non null resultant was found.
+            else:
+                for iRoot in iRootsList:
+                    # A root is given as a (value, multiplicity) tuple.
+                    iRootsSet.add(iRoot[0])
+        # End loop for polyPair in polyParsList. We only loop again if a 
+        # None or zero resultant is found.
         #### Prepare for results for the current interval..
         intervalResultsList = []
         intervalResultsList.append((lb, ub))
         #### Check roots.
         rootsResultsList = []
-        for root in rootsSet:
+        for iRoot in iRootsSet:
             specificRootResultsList = []
-            failingBounds = []
-            intIntPdivN = intIntP(root[0], root[1]) / N
-            if int(intIntPdivN) != intIntPdivN:
-                continue # Next root
+            failingBounds           = []
             # Root qualifies for modular equation, test it for hardness to round.
-            hardToRoundCaseAsFloat = RRR((icAsInt + root[0]) / toIntegerFactor)
+            hardToRoundCaseAsFloat = RRR((icAsInt + iRoot) / toIntegerFactor)
             #print "Before unscaling:", hardToRoundCaseAsFloat.n(prec=precision)
             #print scalingFunction
             scaledHardToRoundCaseAsFloat = \
@@ -1812,26 +2291,18 @@
                            2^-(targetHardnessToRound),
                            RRR):
                 print hardToRoundCaseAsFloat, "is HTRN case."
+                specificRootResultsList.append(hardToRoundCaseAsFloat.n().str(base=2))
                 if lb <= hardToRoundCaseAsFloat and hardToRoundCaseAsFloat <= ub:
                     print "Found in interval."
                 else:
                     print "Found out of interval."
-                specificRootResultsList.append(hardToRoundCaseAsFloat.n().str(base=2))
-                # Check the root is in the bounds
-                if abs(root[0]) > iBound or abs(root[1]) > tBound:
-                    print "Root", root, "is out of bounds for modular equation."
-                    if abs(root[0]) > iBound:
-                        print "root[0]:", root[0]
-                        print "i bound:", iBound
-                        failingBounds.append('i')
-                        failingBounds.append(root[0])
-                        failingBounds.append(iBound)
-                    if abs(root[1]) > tBound:
-                        print "root[1]:", root[1]
-                        print "t bound:", tBound
-                        failingBounds.append('t')
-                        failingBounds.append(root[1])
-                        failingBounds.append(tBound)
+                # Check the i root is within the i bound.
+                if abs(iRoot) > iBound:
+                    print "IRoot", iRoot, "is out of bounds for modular equation."
+                    print "i bound:", iBound
+                    failingBounds.append('i')
+                    failingBounds.append(iRoot)
+                    failingBounds.append(iBound)
                 if len(failingBounds) > 0:
                     specificRootResultsList.append(failingBounds)
             else: # From slz_is_htrn...
@@ -1923,5 +2394,5 @@
     print "Global Wall time:", globalWallTime
     print "Global CPU time:", globalCpuTime
     ## Output counters
-# End srs_compute_lattice_volume
- 
\ No newline at end of file
+# End srs_runSLZ-v04
+
