You can not select more than 25 topics Topics must start with a letter or number, can include dashes ('-') and can be up to 35 characters long.

group_prover.sage 12KB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140141142143144145146147148149150151152153154155156157158159160161162163164165166167168169170171172173174175176177178179180181182183184185186187188189190191192193194195196197198199200201202203204205206207208209210211212213214215216217218219220221222223224225226227228229230231232233234235236237238239240241242243244245246247248249250251252253254255256257258259260261262263264265266267268269270271272273274275276277278279280281282283284285286287288289290291292293294295296297298299300301302303304305306307308309310311312313314315316317318319320321322323324325326327
  1. # This code supports verifying group implementations which have branches
  2. # or conditional statements (like cmovs), by allowing each execution path
  3. # to independently set assumptions on input or intermediary variables.
  4. #
  5. # The general approach is:
  6. # * A constraint is a tuple of two sets of symbolic expressions:
  7. # the first of which are required to evaluate to zero, the second of which
  8. # are required to evaluate to nonzero.
  9. # - A constraint is said to be conflicting if any of its nonzero expressions
  10. # is in the ideal with basis the zero expressions (in other words: when the
  11. # zero expressions imply that one of the nonzero expressions are zero).
  12. # * There is a list of laws that describe the intended behaviour, including
  13. # laws for addition and doubling. Each law is called with the symbolic point
  14. # coordinates as arguments, and returns:
  15. # - A constraint describing the assumptions under which it is applicable,
  16. # called "assumeLaw"
  17. # - A constraint describing the requirements of the law, called "require"
  18. # * Implementations are transliterated into functions that operate as well on
  19. # algebraic input points, and are called once per combination of branches
  20. # executed. Each execution returns:
  21. # - A constraint describing the assumptions this implementation requires
  22. # (such as Z1=1), called "assumeFormula"
  23. # - A constraint describing the assumptions this specific branch requires,
  24. # but which is by construction guaranteed to cover the entire space by
  25. # merging the results from all branches, called "assumeBranch"
  26. # - The result of the computation
  27. # * All combinations of laws with implementation branches are tried, and:
  28. # - If the combination of assumeLaw, assumeFormula, and assumeBranch results
  29. # in a conflict, it means this law does not apply to this branch, and it is
  30. # skipped.
  31. # - For others, we try to prove the require constraints hold, assuming the
  32. # information in assumeLaw + assumeFormula + assumeBranch, and if this does
  33. # not succeed, we fail.
  34. # + To prove an expression is zero, we check whether it belongs to the
  35. # ideal with the assumed zero expressions as basis. This test is exact.
  36. # + To prove an expression is nonzero, we check whether each of its
  37. # factors is contained in the set of nonzero assumptions' factors.
  38. # This test is not exact, so various combinations of original and
  39. # reduced expressions' factors are tried.
  40. # - If we succeed, we print out the assumptions from assumeFormula that
  41. # weren't implied by assumeLaw already. Those from assumeBranch are skipped,
  42. # as we assume that all constraints in it are complementary with each other.
  43. #
  44. # Based on the sage verification scripts used in the Explicit-Formulas Database
  45. # by Tanja Lange and others, see http://hyperelliptic.org/EFD
  46. class fastfrac:
  47. """Fractions over rings."""
  48. def __init__(self,R,top,bot=1):
  49. """Construct a fractional, given a ring, a numerator, and denominator."""
  50. self.R = R
  51. if parent(top) == ZZ or parent(top) == R:
  52. self.top = R(top)
  53. self.bot = R(bot)
  54. elif top.__class__ == fastfrac:
  55. self.top = top.top
  56. self.bot = top.bot * bot
  57. else:
  58. self.top = R(numerator(top))
  59. self.bot = R(denominator(top)) * bot
  60. def iszero(self,I):
  61. """Return whether this fraction is zero given an ideal."""
  62. return self.top in I and self.bot not in I
  63. def reduce(self,assumeZero):
  64. zero = self.R.ideal(list(map(numerator, assumeZero)))
  65. return fastfrac(self.R, zero.reduce(self.top)) / fastfrac(self.R, zero.reduce(self.bot))
  66. def __add__(self,other):
  67. """Add two fractions."""
  68. if parent(other) == ZZ:
  69. return fastfrac(self.R,self.top + self.bot * other,self.bot)
  70. if other.__class__ == fastfrac:
  71. return fastfrac(self.R,self.top * other.bot + self.bot * other.top,self.bot * other.bot)
  72. return NotImplemented
  73. def __sub__(self,other):
  74. """Subtract two fractions."""
  75. if parent(other) == ZZ:
  76. return fastfrac(self.R,self.top - self.bot * other,self.bot)
  77. if other.__class__ == fastfrac:
  78. return fastfrac(self.R,self.top * other.bot - self.bot * other.top,self.bot * other.bot)
  79. return NotImplemented
  80. def __neg__(self):
  81. """Return the negation of a fraction."""
  82. return fastfrac(self.R,-self.top,self.bot)
  83. def __mul__(self,other):
  84. """Multiply two fractions."""
  85. if parent(other) == ZZ:
  86. return fastfrac(self.R,self.top * other,self.bot)
  87. if other.__class__ == fastfrac:
  88. return fastfrac(self.R,self.top * other.top,self.bot * other.bot)
  89. return NotImplemented
  90. def __rmul__(self,other):
  91. """Multiply something else with a fraction."""
  92. return self.__mul__(other)
  93. def __truediv__(self,other):
  94. """Divide two fractions."""
  95. if parent(other) == ZZ:
  96. return fastfrac(self.R,self.top,self.bot * other)
  97. if other.__class__ == fastfrac:
  98. return fastfrac(self.R,self.top * other.bot,self.bot * other.top)
  99. return NotImplemented
  100. # Compatibility wrapper for Sage versions based on Python 2
  101. def __div__(self,other):
  102. """Divide two fractions."""
  103. return self.__truediv__(other)
  104. def __pow__(self,other):
  105. """Compute a power of a fraction."""
  106. if parent(other) == ZZ:
  107. if other < 0:
  108. # Negative powers require flipping top and bottom
  109. return fastfrac(self.R,self.bot ^ (-other),self.top ^ (-other))
  110. else:
  111. return fastfrac(self.R,self.top ^ other,self.bot ^ other)
  112. return NotImplemented
  113. def __str__(self):
  114. return "fastfrac((" + str(self.top) + ") / (" + str(self.bot) + "))"
  115. def __repr__(self):
  116. return "%s" % self
  117. def numerator(self):
  118. return self.top
  119. class constraints:
  120. """A set of constraints, consisting of zero and nonzero expressions.
  121. Constraints can either be used to express knowledge or a requirement.
  122. Both the fields zero and nonzero are maps from expressions to description
  123. strings. The expressions that are the keys in zero are required to be zero,
  124. and the expressions that are the keys in nonzero are required to be nonzero.
  125. Note that (a != 0) and (b != 0) is the same as (a*b != 0), so all keys in
  126. nonzero could be multiplied into a single key. This is often much less
  127. efficient to work with though, so we keep them separate inside the
  128. constraints. This allows higher-level code to do fast checks on the individual
  129. nonzero elements, or combine them if needed for stronger checks.
  130. We can't multiply the different zero elements, as it would suffice for one of
  131. the factors to be zero, instead of all of them. Instead, the zero elements are
  132. typically combined into an ideal first.
  133. """
  134. def __init__(self, **kwargs):
  135. if 'zero' in kwargs:
  136. self.zero = dict(kwargs['zero'])
  137. else:
  138. self.zero = dict()
  139. if 'nonzero' in kwargs:
  140. self.nonzero = dict(kwargs['nonzero'])
  141. else:
  142. self.nonzero = dict()
  143. def negate(self):
  144. return constraints(zero=self.nonzero, nonzero=self.zero)
  145. def __add__(self, other):
  146. zero = self.zero.copy()
  147. zero.update(other.zero)
  148. nonzero = self.nonzero.copy()
  149. nonzero.update(other.nonzero)
  150. return constraints(zero=zero, nonzero=nonzero)
  151. def __str__(self):
  152. return "constraints(zero=%s,nonzero=%s)" % (self.zero, self.nonzero)
  153. def __repr__(self):
  154. return "%s" % self
  155. def conflicts(R, con):
  156. """Check whether any of the passed non-zero assumptions is implied by the zero assumptions"""
  157. zero = R.ideal(list(map(numerator, con.zero)))
  158. if 1 in zero:
  159. return True
  160. # First a cheap check whether any of the individual nonzero terms conflict on
  161. # their own.
  162. for nonzero in con.nonzero:
  163. if nonzero.iszero(zero):
  164. return True
  165. # It can be the case that entries in the nonzero set do not individually
  166. # conflict with the zero set, but their combination does. For example, knowing
  167. # that either x or y is zero is equivalent to having x*y in the zero set.
  168. # Having x or y individually in the nonzero set is not a conflict, but both
  169. # simultaneously is, so that is the right thing to check for.
  170. if reduce(lambda a,b: a * b, con.nonzero, fastfrac(R, 1)).iszero(zero):
  171. return True
  172. return False
  173. def get_nonzero_set(R, assume):
  174. """Calculate a simple set of nonzero expressions"""
  175. zero = R.ideal(list(map(numerator, assume.zero)))
  176. nonzero = set()
  177. for nz in map(numerator, assume.nonzero):
  178. for (f,n) in nz.factor():
  179. nonzero.add(f)
  180. rnz = zero.reduce(nz)
  181. for (f,n) in rnz.factor():
  182. nonzero.add(f)
  183. return nonzero
  184. def prove_nonzero(R, exprs, assume):
  185. """Check whether an expression is provably nonzero, given assumptions"""
  186. zero = R.ideal(list(map(numerator, assume.zero)))
  187. nonzero = get_nonzero_set(R, assume)
  188. expl = set()
  189. ok = True
  190. for expr in exprs:
  191. if numerator(expr) in zero:
  192. return (False, [exprs[expr]])
  193. allexprs = reduce(lambda a,b: numerator(a)*numerator(b), exprs, 1)
  194. for (f, n) in allexprs.factor():
  195. if f not in nonzero:
  196. ok = False
  197. if ok:
  198. return (True, None)
  199. ok = True
  200. for (f, n) in zero.reduce(numerator(allexprs)).factor():
  201. if f not in nonzero:
  202. ok = False
  203. if ok:
  204. return (True, None)
  205. ok = True
  206. for expr in exprs:
  207. for (f,n) in numerator(expr).factor():
  208. if f not in nonzero:
  209. ok = False
  210. if ok:
  211. return (True, None)
  212. ok = True
  213. for expr in exprs:
  214. for (f,n) in zero.reduce(numerator(expr)).factor():
  215. if f not in nonzero:
  216. expl.add(exprs[expr])
  217. if expl:
  218. return (False, list(expl))
  219. else:
  220. return (True, None)
  221. def prove_zero(R, exprs, assume):
  222. """Check whether all of the passed expressions are provably zero, given assumptions"""
  223. r, e = prove_nonzero(R, dict(map(lambda x: (fastfrac(R, x.bot, 1), exprs[x]), exprs)), assume)
  224. if not r:
  225. return (False, map(lambda x: "Possibly zero denominator: %s" % x, e))
  226. zero = R.ideal(list(map(numerator, assume.zero)))
  227. nonzero = prod(x for x in assume.nonzero)
  228. expl = []
  229. for expr in exprs:
  230. if not expr.iszero(zero):
  231. expl.append(exprs[expr])
  232. if not expl:
  233. return (True, None)
  234. return (False, expl)
  235. def describe_extra(R, assume, assumeExtra):
  236. """Describe what assumptions are added, given existing assumptions"""
  237. zerox = assume.zero.copy()
  238. zerox.update(assumeExtra.zero)
  239. zero = R.ideal(list(map(numerator, assume.zero)))
  240. zeroextra = R.ideal(list(map(numerator, zerox)))
  241. nonzero = get_nonzero_set(R, assume)
  242. ret = set()
  243. # Iterate over the extra zero expressions
  244. for base in assumeExtra.zero:
  245. if base not in zero:
  246. add = []
  247. for (f, n) in numerator(base).factor():
  248. if f not in nonzero:
  249. add += ["%s" % f]
  250. if add:
  251. ret.add((" * ".join(add)) + " = 0 [%s]" % assumeExtra.zero[base])
  252. # Iterate over the extra nonzero expressions
  253. for nz in assumeExtra.nonzero:
  254. nzr = zeroextra.reduce(numerator(nz))
  255. if nzr not in zeroextra:
  256. for (f,n) in nzr.factor():
  257. if zeroextra.reduce(f) not in nonzero:
  258. ret.add("%s != 0" % zeroextra.reduce(f))
  259. return ", ".join(x for x in ret)
  260. def check_symbolic(R, assumeLaw, assumeAssert, assumeBranch, require):
  261. """Check a set of zero and nonzero requirements, given a set of zero and nonzero assumptions"""
  262. assume = assumeLaw + assumeAssert + assumeBranch
  263. if conflicts(R, assume):
  264. # This formula does not apply
  265. return None
  266. describe = describe_extra(R, assumeLaw + assumeBranch, assumeAssert)
  267. ok, msg = prove_zero(R, require.zero, assume)
  268. if not ok:
  269. return "FAIL, %s fails (assuming %s)" % (str(msg), describe)
  270. res, expl = prove_nonzero(R, require.nonzero, assume)
  271. if not res:
  272. return "FAIL, %s fails (assuming %s)" % (str(expl), describe)
  273. if describe != "":
  274. return "OK (assuming %s)" % describe
  275. else:
  276. return "OK"
  277. def concrete_verify(c):
  278. for k in c.zero:
  279. if k != 0:
  280. return (False, c.zero[k])
  281. for k in c.nonzero:
  282. if k == 0:
  283. return (False, c.nonzero[k])
  284. return (True, None)