Skip to content

HTTPS clone URL

Subversion checkout URL

You can clone with
or
.
Download ZIP
Browse files

initial drop

  • Loading branch information...
commit 0be4f5c5f0d802660c332b1f1ed9f1cbb64d1fd5 0 parents
@ohpauleez authored
Showing with 34,790 additions and 0 deletions.
  1. +26 −0 COPYRIGHT
  2. +83 −0 RCS/ABS.c,v
  3. +71 −0 RCS/ACLOCK.c,v
  4. +98 −0 RCS/ADDTO.c,v
  5. +47 −0 RCS/AFGEN.c,v
  6. +90 −0 RCS/AFNEG.c,v
  7. +156 −0 RCS/AFPBRI.c,v
  8. +134 −0 RCS/AFPHIIR.c,v
  9. +128 −0 RCS/AFPHIP.c,v
  10. +132 −0 RCS/AFPICR.c,v
  11. +126 −0 RCS/AFPQR.c,v
  12. +101 −0 RCS/AFPRCL.c,v
  13. +306 −0 RCS/AFPRII.c,v
  14. +131 −0 RCS/AFPRLS.c,v
  15. +94 −0 RCS/AFPRRI.c,v
  16. +155 −0 RCS/AFPRRS.c,v
  17. +80 −0 RCS/AFPSIP.c,v
  18. +107 −0 RCS/AFPSIR.c,v
  19. +187 −0 RCS/AFPSUM.c,v
  20. +154 −0 RCS/AFUPBRI.c,v
  21. +133 −0 RCS/AFUPFAC.c,v
  22. +121 −0 RCS/AFUPFMRC.c,v
  23. +100 −0 RCS/AFUPFRPT.c,v
  24. +272 −0 RCS/AFUPGC.c,v
  25. +98 −0 RCS/AFUPGC1.c,v
  26. +425 −0 RCS/AFUPHIBRI.c,v
  27. +91 −0 RCS/AFUPIIR.c,v
  28. +83 −0 RCS/AFUPQR.c,v
  29. +225 −0 RCS/AFUPSFN.c,v
  30. +157 −0 RCS/AFUPSIBRI.c,v
  31. +77 −0 RCS/AGIBL.c,v
  32. +51 −0 RCS/AGICOPY.c,v
  33. +52 −0 RCS/AGIDP2.c,v
  34. +304 −0 RCS/AGIGCD2.c,v
  35. +230 −0 RCS/AGIGCDAE.c,v
  36. +144 −0 RCS/AGIGCDW.c,v
  37. +56 −0 RCS/AGIGI.c,v
  38. +111 −0 RCS/AGIMD.c,v
  39. +51 −0 RCS/AGIMP2.c,v
  40. +105 −0 RCS/AGIMU.c,v
  41. +91 −0 RCS/AGINC.c,v
  42. +93 −0 RCS/AGINORM.c,v
  43. +60 −0 RCS/AGIPROD.c,v
  44. +114 −0 RCS/AGIQHQ.c,v
  45. +88 −0 RCS/AGIRP.c,v
  46. +73 −0 RCS/AGIRSUM.c,v
  47. +52 −0 RCS/AGISUM.c,v
  48. +53 −0 RCS/AGITR.c,v
  49. +53 −0 RCS/AGIWRITE.c,v
  50. +52 −0 RCS/AGIZERO.c,v
  51. +80 −0 RCS/AICOMP.c,v
  52. +57 −0 RCS/AICOPY.c,v
  53. +70 −0 RCS/AIDP2.c,v
  54. +69 −0 RCS/AII.c,v
  55. +83 −0 RCS/AIMP2.c,v
  56. +85 −0 RCS/AINQ.c,v
  57. +133 −0 RCS/AIPROD.c,v
  58. +202 −0 RCS/AIQR.c,v
  59. +162 −0 RCS/AISUM.c,v
  60. +95 −0 RCS/AITR.c,v
  61. +86 −0 RCS/AITRS.c,v
  62. +48 −0 RCS/AIWRITE.c,v
  63. +68 −0 RCS/AIZERO.c,v
  64. +87 −0 RCS/ALSIL.c,v
  65. +70 −0 RCS/AMLM.c,v
  66. +89 −0 RCS/AMSIGN.c,v
  67. +243 −0 RCS/AMSIGNIR.c,v
  68. +75 −0 RCS/AMUPBHT.c,v
  69. +152 −0 RCS/AMUPRINCS.c,v
  70. +115 −0 RCS/ANHI.c,v
  71. +81 −0 RCS/ANR.c,v
  72. +86 −0 RCS/ANSI.c,v
  73. +169 −0 RCS/ARGSACLIB.c,v
  74. +66 −0 RCS/ARIE.c,v
  75. +296 −0 RCS/ASSPR.c,v
  76. +54 −0 RCS/AWCOPY.c,v
  77. +410 −0 RCS/BEGINSACLIB.c,v
  78. +82 −0 RCS/BERNOULLINUM.c,v
  79. +109 −0 RCS/BERNOULLIPOL.c,v
  80. +54 −0 RCS/BRILBRI.c,v
  81. +68 −0 RCS/CHEBY.c,v
  82. +71 −0 RCS/CLEAR.c,v
  83. +95 −0 RCS/COMP2.c,v
  84. +94 −0 RCS/COPYARRAY.c,v
  85. +71 −0 RCS/COPYTO.c,v
  86. +71 −0 RCS/CSFAM.c,v
  87. +55 −0 RCS/CSFS.c,v
  88. +49 −0 RCS/CTMI.c,v
  89. +148 −0 RCS/DDPCC.c,v
  90. +156 −0 RCS/DDRPCC.c,v
  91. +222 −0 RCS/DEGCD.c,v
  92. +250 −0 RCS/DGCD.c,v
  93. +87 −0 RCS/DLINIT.c,v
  94. +68 −0 RCS/DLINV.c,v
  95. +73 −0 RCS/DLNEG.c,v
  96. +108 −0 RCS/DLOG2.c,v
  97. +54 −0 RCS/DLPROD.c,v
  98. +59 −0 RCS/DLSUM.c,v
  99. +175 −0 RCS/DPCC.c,v
  100. +135 −0 RCS/DPGEN.c,v
  101. +165 −0 RCS/DQR.c,v
  102. +122 −0 RCS/DRPCC.c,v
  103. +152 −0 RCS/DSMC.c,v
  104. +85 −0 RCS/EVEN.c,v
  105. +72 −0 RCS/EXPF.c,v
  106. +238 −0 RCS/FAIL.c,v
  107. +86 −0 RCS/FCOMP.c,v
  108. +53 −0 RCS/FCOPY.c,v
  109. +55 −0 RCS/FDIF.c,v
  110. +55 −0 RCS/FIRST5.c,v
  111. +56 −0 RCS/FIRST6.c,v
  112. +57 −0 RCS/FIRST7.c,v
  113. +57 −0 RCS/FIRST8.c,v
  114. +79 −0 RCS/FLBRN.c,v
  115. +106 −0 RCS/FPCATCH.c,v
  116. +79 −0 RCS/FPHAND.c,v
  117. +190 −0 RCS/FPROD.c,v
  118. +97 −0 RCS/FPROD1.c,v
  119. +122 −0 RCS/FPROD2.c,v
  120. +116 −0 RCS/FPROD21.c,v
  121. +133 −0 RCS/FQUOT.c,v
  122. +58 −0 RCS/FRAPCR.c,v
  123. +54 −0 RCS/FRAPFMD.c,v
  124. +58 −0 RCS/FRAPFREE.c,v
  125. +88 −0 RCS/FRAPGET.c,v
  126. +81 −0 RCS/FRAPMON.c,v
  127. +95 −0 RCS/FRAPREM.c,v
  128. +66 −0 RCS/FREEARRAY.c,v
  129. +52 −0 RCS/FREEMATRIX.c,v
  130. +73 −0 RCS/FREINV.c,v
  131. +63 −0 RCS/FREPROD.c,v
  132. +72 −0 RCS/FRUPCR.c,v
  133. +81 −0 RCS/FRUPGCD.c,v
  134. +114 −0 RCS/FSUM.c,v
  135. +195 −0 RCS/FSUMDEOS1.c,v
  136. +300 −0 RCS/FSUMDEOS2.c,v
  137. +509 −0 RCS/FSUMDESS.c,v
  138. +300 −0 RCS/FSUMSEOS.c,v
  139. +120 −0 RCS/FSUMSESS.c,v
  140. +84 −0 RCS/GCAMALLOC.c,v
  141. +57 −0 RCS/GCATL.c,v
  142. +243 −0 RCS/GCSI.c,v
  143. +133 −0 RCS/GDPGEN.c,v
  144. +52 −0 RCS/GETAHIA.c,v
  145. +129 −0 RCS/GETARRAY.c,v
  146. +70 −0 RCS/GETHIPARRAY.c,v
  147. +99 −0 RCS/GETLIST.c,v
  148. +75 −0 RCS/GETMATRIX.c,v
  149. +61 −0 RCS/GIAGI.c,v
  150. +98 −0 RCS/GIBGCD.c,v
  151. +54 −0 RCS/GICONJ.c,v
  152. +66 −0 RCS/GIDIF.c,v
  153. +53 −0 RCS/GIFP.c,v
  154. +62 −0 RCS/GIFQA.c,v
  155. +67 −0 RCS/GIGCD.c,v
  156. +55 −0 RCS/GIHQ.c,v
  157. +86 −0 RCS/GIMS.c,v
  158. +54 −0 RCS/GINEG.c,v
  159. +58 −0 RCS/GINORM.c,v
  160. +58 −0 RCS/GINQ.c,v
  161. +84 −0 RCS/GINQR.c,v
  162. +107 −0 RCS/GIPGEN.c,v
  163. +60 −0 RCS/GIPROD.c,v
  164. +67 −0 RCS/GIRP.c,v
  165. +66 −0 RCS/GISUM.c,v
  166. +73 −0 RCS/GIWRITE.c,v
  167. +87 −0 RCS/GPFUIP.c,v
  168. +98 −0 RCS/GPTOIP.c,v
  169. +110 −0 RCS/GREAD.c,v
  170. +117 −0 RCS/GWRITE.c,v
  171. +52 −0 RCS/HEXP.c,v
  172. +85 −0 RCS/HIACC.c,v
  173. +56 −0 RCS/HIDIF.c,v
  174. +71 −0 RCS/HIDWRITE.c,v
  175. +52 −0 RCS/HILBRI.c,v
  176. +102 −0 RCS/HIPBHT.c,v
  177. +65 −0 RCS/HIPCHT.c,v
  178. +57 −0 RCS/HIPCOPY.c,v
  179. +87 −0 RCS/HIPDWRITE.c,v
  180. +127 −0 RCS/HIPFES.c,v
  181. +61 −0 RCS/HIPIEVAL.c,v
  182. +119 −0 RCS/HIPIR.c,v
  183. +58 −0 RCS/HIPLWRITE.c,v
  184. +62 −0 RCS/HIPNEGT.c,v
  185. +97 −0 RCS/HIPPRB.c,v
  186. +124 −0 RCS/HIPROD.c,v
  187. +225 −0 RCS/HIPRRID.c,v
  188. +283 −0 RCS/HIPRRISD.c,v
  189. +57 −0 RCS/HIPRT.c,v
  190. +88 −0 RCS/HIPSV.c,v
  191. +68 −0 RCS/HIPTR1.c,v
  192. +150 −0 RCS/HIPVCHT.c,v
  193. +98 −0 RCS/HIQUOT.c,v
  194. +57 −0 RCS/HISIGN.c,v
  195. +56 −0 RCS/HISUM.c,v
  196. +55 −0 RCS/HSIGN.c,v
  197. +66 −0 RCS/IACOMPA.c,v
  198. +96 −0 RCS/IAI.c,v
  199. +84 −0 RCS/IBCIND.c,v
  200. +64 −0 RCS/IBCOEFS.c,v
  201. +69 −0 RCS/IBPPOL.c,v
  202. +156 −0 RCS/IBPPOS.c,v
  203. +70 −0 RCS/ICRAND.c,v
  204. +111 −0 RCS/IDEGCD.c,v
  205. +56 −0 RCS/IDENTMAT.c,v
  206. +188 −0 RCS/IDEQ.c,v
  207. +300 −0 RCS/IDIF.c,v
  208. +339 −0 RCS/IDIFA.c,v
  209. +156 −0 RCS/IDLCOMB.c,v
  210. +171 −0 RCS/IDP2.c,v
  211. +203 −0 RCS/IDPR.c,v
  212. +71 −0 RCS/IDQ.c,v
  213. +219 −0 RCS/IDQR.c,v
  214. +120 −0 RCS/IDQRA.c,v
  215. +71 −0 RCS/IDREM.c,v
  216. +58 −0 RCS/IEAS.c,v
  217. +55 −0 RCS/IEEEDWRITE.c,v
  218. +90 −0 RCS/IEEELBRN.c,v
  219. +124 −0 RCS/IEEENEIGH.c,v
  220. +99 −0 RCS/IEEEROUND.c,v
  221. +69 −0 RCS/IEEEWRITE.c,v
  222. +99 −0 RCS/IEGCD.c,v
  223. +341 −0 RCS/IEQ.c,v
  224. +137 −0 RCS/IF.c,v
  225. +256 −0 RCS/IFATL.c,v
  226. +157 −0 RCS/IFCL2.c,v
  227. +100 −0 RCS/IFEC.c,v
  228. +336 −0 RCS/IFLTA.c,v
  229. +161 −0 RCS/IGCD.c,v
  230. +98 −0 RCS/IGCDCF.c,v
  231. +63 −0 RCS/IHDREM.c,v
  232. +79 −0 RCS/IHEAS.c,v
  233. +189 −0 RCS/IHEGCD.c,v
  234. +82 −0 RCS/IHI.c,v
  235. +212 −0 RCS/IIEEE.c,v
  236. +163 −0 RCS/IIEEET.c,v
  237. +52 −0 RCS/ILBRN.c,v
  238. +99 −0 RCS/ILCM.c,v
  239. +38 −0 RCS/ILENGTH.c,v
  240. +41 −0 RCS/ILOG2A.c,v
  241. +59 −0 RCS/ILOGB.c,v
  242. +58 −0 RCS/ILSA.c,v
  243. +109 −0 RCS/ILWORD.c,v
  244. +101 −0 RCS/ILWORDS.c,v
  245. +255 −0 RCS/IMP2.c,v
  246. +146 −0 RCS/IMPB.c,v
  247. +217 −0 RCS/IMPBAA.c,v
  248. +116 −0 RCS/INEG.c,v
  249. +51 −0 RCS/INEGA.c,v
  250. +96 −0 RCS/INFOSACLIB.c,v
  251. +83 −0 RCS/INQ.c,v
  252. +117 −0 RCS/INQR.c,v
  253. +59 −0 RCS/INSET.c,v
  254. +55 −0 RCS/INVPERM.c,v
  255. +122 −0 RCS/IORD2.c,v
  256. +93 −0 RCS/IPBEZM.c,v
  257. +122 −0 RCS/IPBHT.c,v
  258. +127 −0 RCS/IPBHTLV.c,v
  259. +141 −0 RCS/IPBHTMV.c,v
  260. +112 −0 RCS/IPCA.c,v
  261. +166 −0 RCS/IPCEVP.c,v
  262. +90 −0 RCS/IPDQNR.c,v
  263. +90 −0 RCS/IPDSCR.c,v
  264. +58 −0 RCS/IPDSCRBEZ.c,v
  265. +218 −0 RCS/IPEQ.c,v
  266. +86 −0 RCS/IPEXP.c,v
  267. +175 −0 RCS/IPFAC.c,v
  268. +246 −0 RCS/IPFLC.c,v
  269. +435 −0 RCS/IPGCDC.c,v
  270. +87 −0 RCS/IPGFCB.c,v
  271. +111 −0 RCS/IPGSUB.c,v
  272. +76 −0 RCS/IPGTDRAN.c,v
  273. +78 −0 RCS/IPIIRB.c,v
  274. +304 −0 RCS/IPIIS.c,v
  275. +73 −0 RCS/IPIISS.c,v
  276. +114 −0 RCS/IPINT.c,v
  277. +128 −0 RCS/IPIP.c,v
  278. +113 −0 RCS/IPIQ.c,v
  279. +191 −0 RCS/IPIQH.c,v
  280. +58 −0 RCS/IPLEV.c,v
  281. +76 −0 RCS/IPLSEVAL.c,v
  282. +93 −0 RCS/IPLSILM.c,v
  283. +96 −0 RCS/IPMONFACT.c,v
  284. +197 −0 RCS/IPP2P.c,v
  285. +81 −0 RCS/IPPFAC2.c,v
  286. +72 −0 RCS/IPPNPRS.c,v
  287. +208 −0 RCS/IPPSC.c,v
  288. +97 −0 RCS/IPPVED.c,v
  289. +86 −0 RCS/IPQ.c,v
  290. +82 −0 RCS/IPRAN.c,v
  291. +124 −0 RCS/IPRCH.c,v
  292. +132 −0 RCS/IPRCHS.c,v
  293. +258 −0 RCS/IPRES.c,v
  294. +97 −0 RCS/IPRESBEZ.c,v
  295. +135 −0 RCS/IPRESPRS.c,v
  296. +80 −0 RCS/IPRIP.c,v
  297. +61 −0 RCS/IPRIST.c,v
  298. +55 −0 RCS/IPRNME.c,v
  299. +541 −0 RCS/IPROD.c,v
  300. +120 −0 RCS/IPROD2X2.c,v
Sorry, we could not display the entire diff because too many files (2,226) changed.
26 COPYRIGHT
@@ -0,0 +1,26 @@
+ SACLIB (C) Copyright 1993 by Kurt Goedel Institute
+
+
+
+The SACLIB system source code and User's Guide are made available
+free of charge by the Kurt Goedel Institute on behalf of the SACLIB
+Group.
+
+Persons or institutions receiving it are pledged not to distribute it to
+others. Instead, individuals wishing to acquire the system should obtain it
+by ftp directly from the Kurt Goedel Institute, informing the Institute of
+the acquisition. Thereby the SACLIB Group will know who has the system and
+be able to inform all users of any corrections or newer versions.
+
+Users are kindly asked to cite their use of the system in any resulting
+publications or in any application packages built upon SACLIB.
+
+Neither SACLIB nor any part thereof may be incorporated in any commercial
+software product without the consent of the authors. Users developing
+non-commercial application packages are kindly asked to inform us.
+
+Requests or proposals for changes or additions to the system will
+be welcomed and given consideration.
+
+SACLIB is offered without warranty of any kind, either expressed or implied.
+However reports of bugs or problems are encouraged.
83 RCS/ABS.c,v
@@ -0,0 +1,83 @@
+head 1.2;
+access;
+symbols;
+locks
+ saclib:1.2; strict;
+comment @ * @;
+
+
+1.2
+date 94.08.21.11.56.29; author George; state spec;
+branches;
+next 1.1;
+
+1.1
+date 94.08.21.11.55.37; author saclib; state init;
+branches;
+next ;
+
+
+desc
+@Absolute value.
+@
+
+
+1.2
+log
+@GAMMA integer was replaced by C integer.
+@
+text
+@/*===========================================================================
+ b <- ABS(a)
+
+Absolute value.
+
+Inputs
+ a : a C integer.
+
+Outputs
+ b : the absolute value of a.
+===========================================================================*/
+#include "saclib.h"
+
+Word ABS(a)
+ Word a;
+{
+ Word b;
+
+Step1: /* Get the absolute value. */
+ if (a >= 0)
+ b = a;
+ else
+ b = -a;
+
+Return: /* Prepare for return. */
+ return(b);
+}
+@
+
+
+1.1
+log
+@Initial revision
+@
+text
+@d1 2
+a2 2
+/*======================================================================
+ b <- ABS(a)
+d7 1
+a7 1
+ a : a gamma-integer.
+d11 1
+a11 1
+======================================================================*/
+a17 1
+ /* hide algorithm */
+d21 1
+a21 1
+ b = a;
+d23 1
+a23 1
+ b = -a;
+@
71 RCS/ACLOCK.c,v
@@ -0,0 +1,71 @@
+head 1.2;
+access;
+symbols;
+locks
+ saclib:1.2; strict;
+comment @ * @;
+
+
+1.2
+date 94.08.21.11.57.50; author George; state spec;
+branches;
+next 1.1;
+
+1.1
+date 94.08.21.11.57.20; author saclib; state init;
+branches;
+next ;
+
+
+desc
+@Clock minus garbage collection time.
+@
+
+
+1.2
+log
+@GAMMA integer was replaced by C integer.
+@
+text
+@/*===========================================================================
+ T <- ACLOCK()
+
+Clock minus garbage collection time.
+
+Outputs
+ T : a C integer, the system clock time minus garbage collection
+ time in milliseconds.
+===========================================================================*/
+#include "saclib.h"
+
+Word ACLOCK()
+{
+ Word T;
+
+Step1: /* Compute. */
+ T = CLOCK() - TAU;
+
+Return: /* Prepare for return. */
+ return(T);
+}
+@
+
+
+1.1
+log
+@Initial revision
+@
+text
+@d1 2
+a2 2
+/*======================================================================
+ T <- ACLOCK()
+d7 1
+a7 1
+ T : a gamma integer, the system clock time minus garbage collection
+d9 1
+a9 1
+======================================================================*/
+a14 1
+ /* hide algorithm */
+@
98 RCS/ADDTO.c,v
@@ -0,0 +1,98 @@
+head 1.2;
+access;
+symbols;
+locks
+ saclib:1.2; strict;
+comment @ * @;
+
+
+1.2
+date 95.08.30.11.18.31; author Mark; state bug;
+branches;
+next 1.1;
+
+1.1
+date 93.12.13.20.08.39; author Werner; state new;
+branches;
+next ;
+
+
+desc
+@Add array to array.
+@
+
+
+1.2
+log
+@The BDigit declarations were changed to Word since BDigit caused problems
+with isac.
+@
+text
+@/*===========================================================================
+ ADDTO(A,B,n)
+
+Add array to array.
+
+Inputs
+ A, B : arrays of length >= n containing non-negative BETA-digits.
+ n : positive BETA-digit.
+
+Effect : B[0],...,B[n-1] is added to A[0],...,A[n-1],... (result in A).
+===========================================================================*/
+#include "saclib.h"
+
+void ADDTO(A,B,n)
+ Word *A,*B,n;
+{
+ Word c,cp,i;
+
+Step1: /* Add. */
+ cp = 0;
+ for (i = 0; i < n; i++)
+ {
+ c = A[i] + B[i] + cp;
+ if (c >= BETA)
+ {
+ c = c - BETA;
+ cp = 1;
+ }
+ else
+ cp = 0;
+ A[i] = c;
+ }
+
+Step2: /* Propagate carry. */
+ while (cp != 0)
+ {
+ c = A[i] + cp;
+ if (c >= BETA)
+ {
+ c = c - BETA;
+ cp = 1;
+ }
+ else
+ cp =0;
+ A[i] = c;
+ i++;
+ }
+
+Return: /* Prepare for return. */
+ return;
+}
+@
+
+
+1.1
+log
+@Initial revision
+@
+text
+@d15 1
+a15 2
+ Word *A,*B;
+ BDigit n;
+d17 1
+a17 2
+ BDigit cp,i;
+ Word c;
+@
47 RCS/AFGEN.c,v
@@ -0,0 +1,47 @@
+head 1.1;
+access;
+symbols;
+locks
+ saclib:1.1; strict;
+comment @ * @;
+
+
+1.1
+date 94.07.02.23.19.57; author saclib; state new;
+branches;
+next ;
+
+
+desc
+@Algebraic number field generator.
+@
+
+
+1.1
+log
+@Initial revision
+@
+text
+@/*===========================================================================
+ a <- AFGEN()
+
+Algebraic number field generator.
+
+Output
+ a : in Q(alpha). If alpha is any irrational algebraic number, a is
+ the representation of alpha as an element of the algebraic number
+ field Q(alpha).
+===========================================================================*/
+#include "saclib.h"
+
+Word AFGEN()
+{
+ Word a;
+
+Step1: /* Construct. */
+ a = LIST2(RNINT(1),PMON(1,1));
+
+Return: /* Prepare for return. */
+ return(a);
+}
+@
90 RCS/AFNEG.c,v
@@ -0,0 +1,90 @@
+head 1.2;
+access;
+symbols;
+locks
+ saclib:1.2; strict;
+comment @ * @;
+
+
+1.2
+date 94.05.11.12.40.04; author Mark; state typo;
+branches;
+next 1.1;
+
+1.1
+date 94.05.11.12.38.03; author saclib; state init;
+branches;
+next ;
+
+
+desc
+@Algebraic number field negative.
+@
+
+
+1.2
+log
+@Input was labelled Output and vice versa.
+@
+text
+@/*===========================================================================
+ b <- AFNEG(a)
+
+Algebraic number field negation.
+
+Input
+ a : in Q(alpha).
+
+Output
+ b : in Q(alpha), b = -a.
+===========================================================================*/
+#include "saclib.h"
+
+Word AFNEG(a)
+ Word a;
+{
+ Word A,b,r;
+
+Step1: /* Change sign. */
+ if (a == 0) {
+ b = 0;
+ goto Return; }
+ FIRST2(a,&r,&A);
+ b = LIST2(RNNEG(r),A);
+
+Return: /* Prepare for return. */
+ return(b);
+}
+@
+
+
+1.1
+log
+@Initial revision
+@
+text
+@d1 2
+a2 2
+/*======================================================================
+ b <- AFNEG(a)
+d4 1
+a4 1
+Algebraic number field negative.
+d6 1
+a6 1
+Outputs
+d9 1
+a9 1
+Input:
+d11 1
+a11 1
+======================================================================*/
+d19 4
+a22 6
+Step1:
+ if (a == 0)
+ {
+ b = 0;
+ goto Return;
+ }
+@
156 RCS/AFPBRI.c,v
@@ -0,0 +1,156 @@
+head 1.2;
+access;
+symbols;
+locks
+ saclib:1.2; strict;
+comment @ * @;
+
+
+1.2
+date 95.08.09.13.49.10; author Mark; state del;
+branches;
+next 1.1;
+
+1.1
+date 95.08.09.13.42.52; author saclib; state init;
+branches;
+next ;
+
+
+desc
+@Algebraic number field polynomial basis real root isolation.
+@
+
+
+1.2
+log
+@Replaced by AFUPBRI.
+@
+text
+@[ removed from library ]
+/*======================================================================
+ N <- AFPBRI(M,I,L)
+
+Algebraic number field polynomial basis real root isolation.
+
+Inputs
+ M : in Z[x], is the minimal polynomial of alpha.
+ I : an acceptable isolating interval for alpha.
+ L : a list (A1,...,An) of elements of Q(alpha)[X].
+ L is a nonempty squarefree basis.
+
+Outputs
+ N : a list (I1,B1,...,Im,Bm), m >= 0, where I1 < I2 < ... < Im
+ are strongly disjoint isolating intervals for all the real roots of
+ A = A1 * ... * An. Each Ii has binary rational endpoints and is
+ left open and right closed. Bi is the unique Aj which has a root
+ in Ii.
+======================================================================*/
+#include "saclib.h"
+
+Word AFPBRI(M,I,L)
+ Word M,I,L;
+{
+ Word A1,A2,B1,I1,I11,I21,Lp,Lpp,N,S,S1,S2,Sp,Spp,Ss1,Ss2,T,T1,
+ T2,m,s;
+ /* hide A1,A2,B1,Lp,Lpp,Sp,Spp,m,s; */
+
+Step1: /* Isolate roots of each Ai. */
+ Lp = L;
+ S = NIL;
+ do
+ {
+ ADV(Lp,&A1,&Lp);
+ S1 = AFPRCL(M,I,A1);
+ S = COMP(S1,S);
+ }
+ while (Lp != NIL);
+ S = INV(S);
+
+Step2: /* Refine to disjoint isolating intervals. */
+ Lp = L;
+ Sp = S;
+ while (RED(Lp) != NIL)
+ {
+ A1 = FIRST(Lp);
+ S1 = FIRST(Sp);
+ Lpp = RED(Lp);
+ Spp = RED(Sp);
+ do
+ {
+ A2 = FIRST(Lpp);
+ S2 = FIRST(Spp);
+ AFPRLS(M,I,A1,A2,S1,S2,&Ss1,&Ss2);
+ S1 = Ss1;
+ SFIRST(Spp,Ss2);
+ Lpp = RED(Lpp);
+ Spp = RED(Spp);
+ }
+ while (Lpp != NIL);
+ SFIRST(Sp,S1);
+ Lp = RED(Lp);
+ Sp = RED(Sp);
+ }
+
+Step3: /* Prepare to merge intervals. */
+ Lp = L;
+ Sp = S;
+ T = NIL;
+ do
+ {
+ ADV(Lp,&A1,&Lp);
+ ADV(Sp,&S1,&Sp);
+ T1 = NIL;
+ while (S1 != NIL)
+ {
+ ADV2(S1,&I11,&m,&S1);
+ T1 = COMP2(A1,I11,T1);
+ }
+ T = COMP(INV(T1),T);
+ }
+ while (Lp != NIL);
+ T = INV(T);
+
+Step4: /* Merge-sort isolating intervals. */
+ while (RED(T) != NIL)
+ {
+ S = NIL;
+ while (T != NIL && RED(T) != NIL)
+ {
+ ADV2(T,&T1,&T2,&T);
+ S1 = NIL;
+ while (T1 != NIL && T2 != NIL)
+ {
+ I11 = FIRST(T1);
+ I21 = FIRST(T2);
+ s = RNCOMP(FIRST(I11),FIRST(I21));
+ if (s < 0)
+ ADV2(T1,&I1,&B1,&T1);
+ else
+ ADV2(T2,&I1,&B1,&T2);
+ S1 = COMP2(B1,I1,S1);
+ }
+ if (T1 == NIL)
+ T1 = T2;
+ S1 = CONC(INV(S1),T1);
+ S = COMP(S1,S);
+ }
+ if (T != NIL)
+ S = COMP(FIRST(T),S);
+ T = INV(S);
+ }
+ N = FIRST(T);
+
+Return: /* Prepare for return. */
+ return(N);
+}
+@
+
+
+1.1
+log
+@Initial revision
+@
+text
+@d1 1
+@
134 RCS/AFPHIIR.c,v
@@ -0,0 +1,134 @@
+head 1.2;
+access;
+symbols;
+locks
+ saclib:1.2; strict;
+comment @ * @;
+
+
+1.2
+date 2001.07.02.19.02.04; author George; state bug;
+branches;
+next 1.1;
+
+1.1
+date 98.04.23.20.25.52; author George; state new;
+branches;
+next ;
+
+
+desc
+@*
+Algebraic field polynomial hardware interval isolation and
+refinement.
+@
+
+
+1.2
+log
+@Corrected the call of HIPRRID in Step 2.
+@
+text
+@/*======================================================================
+ L <- AFPHIIR(n,I,A,k)
+
+Algebraic field polynomial hardware interval isolation and
+refinement.
+
+Inputs
+ n : a positive beta-digit.
+ I : a hardware interval containing an algebraic number \alpha
+ of degree n.
+ A : a univariate polynomial of positive degree over the
+ algebraic number field Q(\alpha). Each nonzero coefficient
+ of A is represented by a pair (r,a) where r is a
+ rational number and a is an integral polynomial of
+ degree less than n,
+ k : a beta-digit.
+
+Output
+ L : either 0 or a list of isolating intervals for the real
+ roots of A. A(x) is converted to a hardware interval
+ polynomial B(x) using I. If the real roots of B(x)
+ cannot be isolated L = 0, a failure condition.
+ Otherwise L is a list of logarithmic binary rational
+ intervals for A(x). Each interval is either a one-point
+ interval, an interval of width 1 / 2^k, or wider in
+ case width 1 / 2^k can not be achieved using hardware
+ interval airthmetic.
+
+Remark
+ AFPHIIR is floating-point overflow-underflow protected.
+======================================================================*/
+#include "saclib.h"
+
+Word AFPHIIR(n,I,A,k)
+ BDigit n;
+ interval I;
+ Word A;
+ BDigit k;
+{
+ interval *B,Jp,K;
+ BDigit h,j,m,t,T,u;
+ Word J,Kp,L,Lp,w;
+
+Step1: /* Convert A to an interval polynomial B. */
+ AFPHIP(I,A,&B,&T);
+ if (T == 0) {
+ L = 0;
+ goto Return; }
+
+Step2: /* Isolate the roots of B. */
+ m = PDEG(A);
+ HIPRRID(m,B, &L,&T);
+ if (T == 0) {
+ L = 0;
+ goto Return; }
+
+Step3: /* Determine the trend of the first interval. */
+ if (B[m].left > 0) {
+ if (EVEN(m))
+ t = -1;
+ else
+ t = +1; }
+ else {
+ if (EVEN(m))
+ t = +1;
+ else
+ t = -1; }
+
+Step4: /* Refine the isolating intervals. */
+ Lp = L;
+ while (Lp != NIL) {
+ J = FIRST(Lp);
+ w = LBRIW(J);
+ if (w != 0) {
+ h = SECOND(w);
+ if (h < k) {
+ LBRIHI(J,&Jp,&u);
+ if (u == 0) {
+ L = 0;
+ goto Return; }
+ HIPIR(m,B,Jp,t,h,k,&K,&j);
+ Kp = HILBRI(K);
+ SFIRST(Lp,Kp); } }
+ Lp = RED(Lp);
+ t = - t; }
+
+Return: /* Return L. */
+ return(L);
+}
+@
+
+
+1.1
+log
+@Initial revision
+@
+text
+@d52 4
+a55 3
+ L = HIPRRID(m,B);
+ if (L == 0)
+ goto Return;
+@
128 RCS/AFPHIP.c,v
@@ -0,0 +1,128 @@
+head 1.3;
+access;
+symbols;
+locks
+ saclib:1.3; strict;
+comment @ * @;
+
+
+1.3
+date 2000.04.10.19.26.49; author George; state embe;
+branches;
+next 1.2;
+
+1.2
+date 99.04.21.17.32.45; author George; state bug;
+branches;
+next 1.1;
+
+1.1
+date 98.04.23.20.26.59; author George; state new;
+branches;
+next ;
+
+
+desc
+@*
+Algebraic field polynomial to hardware interval polynomial.
+@
+
+
+1.3
+log
+@Improved the specification.
+@
+text
+@/*======================================================================
+ AFPHIP(I,A; B,t)
+
+Algebraic field polynomial to hardware interval polynomial.
+
+Inputs
+ I : a hardware interval containing an algebraic number \alpha.
+ A : a univariate polynomial of positive degree over the
+ algebraic number field Q(\alpha). Each nonzero coefficient
+ of A is represented by a pair (r,a) where r is a
+ rational number and a is an integral polynomial of
+ degree less than the degree of \alpha.
+
+Outputs
+ B : if t = 1, B is a hardware interval polynomial containing A(x).
+ Otherwise t = 0 and the value of B is undefined due to
+ exponent limitation.
+ t : 0 or 1.
+======================================================================*/
+#include "saclib.h"
+
+void AFPHIP(I,A,B_,t_)
+ interval I;
+ Word A;
+ interval **B_;
+ BDigit *t_;
+
+{
+ interval *B;
+ BDigit i,m,t;
+ Word a,Ap;
+ Word V,v;
+
+Step1: /* Get an array for B. */
+ m = PDEG(A);
+ B = GETHIPARRAY(m);
+
+Step2: /* Convert and evaluate coefficients. */
+ Ap = A;
+ for (i = m; i >= 0; i--) {
+ if (Ap == NIL || PDEG(Ap) < i)
+ IHI(0,&B[i],&t);
+ else {
+ Ap = RED(Ap);
+ ADV(Ap,&a,&Ap);
+ ANHI(I,a,&B[i],&t);
+ if (t == 0)
+ goto Return; } }
+
+Return: /* Return B and t. */
+ *B_ = B;
+ *t_ = t;
+ return;
+}
+@
+
+
+1.2
+log
+@In Step2, the line
+ if (Ap == NIL | PDEG(Ap) < i)
+was changed to
+ if (Ap == NIL || PDEG(Ap) < i)
+@
+text
+@d15 4
+a18 6
+ B : a hardware interval polynomial containing A(x), if t = 1;
+ otherwise undefined.
+ t : t = 1 if AFPHIP succeeds; otherwise t = 0.
+
+Remark
+ AFPHIP is floating-point overflow-underflow protected.
+d47 2
+a48 3
+ if (t == 0) {
+ SWRITE("ANHI failure.\n");
+ goto Return; } } }
+@
+
+
+1.1
+log
+@Initial revision
+@
+text
+@d43 1
+a43 1
+ if (Ap == NIL | PDEG(Ap) < i)
+d50 1
+a50 1
+ SWRITE("ANHHI failure.\n");
+@
132 RCS/AFPICR.c,v
@@ -0,0 +1,132 @@
+head 1.2;
+access;
+symbols;
+locks
+ saclib:1.2; strict;
+comment @ * @;
+
+
+1.2
+date 94.06.29.14.57.13; author saclib; state embe;
+branches;
+next 1.1;
+
+1.1
+date 94.06.29.14.56.43; author saclib; state init;
+branches;
+next ;
+
+
+desc
+@Algebraic number field polynomial inverse convert representation.
+@
+
+
+1.2
+log
+@Made some cosmetic changes.
+@
+text
+@/*===========================================================================
+ B <- AFPICR(r,A)
+
+Algebraic number field polynomial inverse convert representation.
+
+Inputs
+ r : a non-negative BETA-integer.
+ A : an element of Q(alpha)[x_1,...,x_r], where the coefficients of A are
+ represented as (a_i,Ai), where a_i is an element of Q and Ai is an
+ element of Z[x].
+
+Outputs
+ B : an element of Q(alpha)[x_1,...,x_r] represented as an element
+ of Q[x,x_1,...,x_r].
+===========================================================================*/
+#include "saclib.h"
+
+Word AFPICR(r,A)
+
+ Word r,A;
+{
+ Word B;
+ Word Ap,a,b,e,rp;
+
+Step1: /* A = 0. */
+ if (A == 0) {
+ B = 0;
+ goto Return; }
+
+Step2: /* r = 0. */
+ if (r == 0) {
+ B = AFICR(A);
+ goto Return; }
+
+Step3: /* r > 0. */
+ rp = r - 1;
+ Ap = A;
+ B = NIL;
+ while (Ap != NIL) {
+ ADV2(Ap,&e,&a,&Ap);
+ b = AFPICR(rp,a);
+ B = COMP2(b,e,B); }
+ B = INV(B);
+
+Return: /* Prepare for return. */
+ return(B);
+}
+@
+
+
+1.1
+log
+@Initial revision
+@
+text
+@d1 2
+a2 2
+/*======================================================================
+ B <- AFPICR(r,A)
+d8 1
+a8 1
+ A : an element of Q(alpha)[X1,...,Xr], where the coefficients of A are
+d10 1
+a10 1
+ element of Z[X].
+d13 3
+a15 3
+ B : an element of Q(alpha)[X1,...,Xr] represented as an element of
+ Q[X,X1,...,Xr].
+======================================================================*/
+d25 9
+a33 13
+Step1: /* A=0. */
+ if (A == 0)
+ {
+ B = 0;
+ goto Return;
+ }
+
+Step2: /* r = 0? */
+ if (r == 0)
+ {
+ B = AFICR(A);
+ goto Return;
+ }
+d35 1
+a35 1
+Step3:
+d39 4
+a42 9
+ while (Ap != NIL)
+ {
+ ADV2(Ap,&e,&a,&Ap);
+ if (rp == 0)
+ b = AFICR(a);
+ else
+ b = AFPICR(rp,a);
+ B = COMP2(b,e,B);
+ }
+d45 1
+a45 1
+Return:
+@
126 RCS/AFPQR.c,v
@@ -0,0 +1,126 @@
+head 1.2;
+access;
+symbols;
+locks
+ saclib:1.2; strict;
+comment @ * @;
+
+
+1.2
+date 93.05.17.14.59.35; author George; state spec;
+branches;
+next 1.1;
+
+1.1
+date 93.05.17.14.55.15; author saclib; state init;
+branches;
+next ;
+
+
+desc
+@Algebraic number field polynomial quotient and remainder.
+@
+
+
+1.2
+log
+@The specification previously defined Q as B / A.
+@
+text
+@/*==========================================================================
+ AFPQR(r,M,A,B; Q,R)
+
+Algebraic number field polynomial quotient and remainder.
+
+Inputs
+ r : a BETA-digit, r >= 1.
+ M : in Z[x], the minimal polynomial for an algebraic number alpha.
+ A,B : in Q(alpha)[X1,...,Xr], B not zero.
+
+Outputs
+ Q,R : in Q(alpha)[X1,...,Xr].
+ Q and R are the unique algebraic number field polynomials
+ such that either B divides A, Q = A / B, and R = 0 or else B
+ does not divide A and A = BQ+R with degree(R) minimal.
+==========================================================================*/
+#include "saclib.h"
+
+void AFPQR(r,M,A,B, Q_,R_)
+ Word r,M,A,B, *Q_,*R_;
+{
+ Word Bp,Q,Q1,Qp,R,Rp,a,b,d,m,n,q,rp,s;
+
+Step1: /* Initialize. */
+ n = PDEG(B);
+ b = PLDCF(B);
+ Bp = PRED(B);
+ Q = NIL;
+ R = A;
+ rp = r - 1;
+
+Step2: /* Compute quotient terms. */
+ while (R != 0) {
+ m = PDEG(R);
+ d = m - n;
+ if (d < 0)
+ goto Step3;
+ a = PLDCF(R);
+ if (rp == 0) {
+ q = AFQ(M,a,b);
+ s = 0;
+ } else
+ AFPQR(rp,M,a,b,&q,&s);
+ if (s != 0)
+ goto Step3;
+ Q = COMP2(q,d,Q);
+ Q1 = LIST2(d,q);
+ Rp = PRED(R);
+ Qp = AFPPR(r,M,Bp,Q1);
+ R = AFPDIF(r,Rp,Qp);
+ }
+
+Step3: /* Finish. */
+ if (Q == NIL)
+ Q = 0;
+ else
+ Q = INV(Q);
+
+Return: /* Prepare for return. */
+ *Q_ = Q;
+ *R_ = R;
+ return;
+}
+@
+
+
+1.1
+log
+@Initial revision
+@
+text
+@d1 1
+a1 1
+/*======================================================================
+d14 3
+a16 3
+ such that either B divides A, Q = B/A, and R = 0 or else B does not
+ divide A and A = BQ+R with degree(R) minimal.
+======================================================================*/
+a22 1
+ /* hide Rp,a,d,m,n,rp; */
+d33 1
+a33 2
+ while (R != 0)
+ {
+d39 1
+a39 2
+ if (rp == 0)
+ {
+d42 1
+a42 2
+ }
+ else
+d51 1
+a51 1
+ }
+@
101 RCS/AFPRCL.c,v
@@ -0,0 +1,101 @@
+head 1.2;
+access;
+symbols;
+locks
+ saclib:1.2; strict;
+comment @ * @;
+
+
+1.2
+date 95.08.09.13.49.47; author Mark; state del;
+branches;
+next 1.1;
+
+1.1
+date 95.08.09.13.43.14; author saclib; state init;
+branches;
+next ;
+
+
+desc
+@Algebraic number field polynomial real root isolation, collins-loos
+algorithm.
+@
+
+
+1.2
+log
+@Replaced by AFUPRICL.
+@
+text
+@[ removed from library ]
+/*======================================================================
+ L <- AFPRCL(M,I,A)
+
+Algebraic number field polynomial real root isolation, collins-loos
+algorithm.
+
+Inputs
+ M : in Z[x], the minimal polynomial of an algebraic number alpha.
+ I : an acceptable isolating interval for alpha.
+ A : in Q(alpha)[x], A is monic and deg(A)=n.
+
+Outputs
+ L : a strong isolation list for the real roots of A.
+======================================================================*/
+#include "saclib.h"
+
+Word AFPRCL(M,I,A)
+ Word M,I,A;
+{
+ Word Ap,As,C,L,Lp,b,d;
+
+Step1: /* Degree zero. */
+ if (PDEG(A) == 0)
+ {
+ L = NIL;
+ goto Return;
+ }
+
+Step2: /* Initialize. */
+ As = A;
+ Ap = AFPDMV(1,M,A);
+ C = NIL;
+
+Step3: /* Compute derivatives. */
+ while (PDEG(As) > 1)
+ {
+ b = PTBCF(1,As);
+ C = COMP(b,C);
+ As = Ap;
+ Ap = AFPDMV(1,M,As);
+ }
+
+Step4: /* Isolate roots. */
+ Lp = NIL;
+ d = AFUPRB(M,I,A);
+ do
+ {
+ L = AFPRII(M,I,As,Ap,d,Lp);
+ if (C == NIL)
+ goto Return;
+ ADV(C,&b,&C);
+ Ap = As;
+ As = AFPINT(1,M,Ap,b);
+ Lp = L;
+ }
+ while (1);
+
+Return: /* Prepare for return. */
+ return(L);
+}
+@
+
+
+1.1
+log
+@Initial revision
+@
+text
+@d1 1
+@
306 RCS/AFPRII.c,v
@@ -0,0 +1,306 @@
+head 1.2;
+access;
+symbols;
+locks
+ saclib:1.2; strict;
+comment @ * @;
+
+
+1.2
+date 95.08.09.13.50.36; author Mark; state del;
+branches;
+next 1.1;
+
+1.1
+date 95.08.09.13.43.45; author saclib; state init;
+branches;
+next ;
+
+
+desc
+@Algebraic number field polynomial real root isolation induction.
+@
+
+
+1.2
+log
+@Replaced by AFUPRII.
+@
+text
+@[ removed from library ]
+/*======================================================================
+ L <- AFPRII(M,J,A,Ap,d,Lp)
+
+Algebraic number field polynomial real root isolation induction.
+
+Inputs
+ M : in Z[x], the minimal polynomial of an algebraic number alpha.
+ J : an acceptable isolating interval for alpha.
+ A : in Q(alpha)[x], positive and deg(A) > 0.
+ Ap : in Q(alpha)[x], the derivative of A.
+ d : a binary rational real root bound for A.
+ Lp : a strong isolation list for Ap.
+
+Outputs
+ L : a strong isolation list for A.
+======================================================================*/
+#include "saclib.h"
+
+Word AFPRII(M,J,A,Ap,d,Lp)
+ Word M,J,A,Ap,d,Lp;
+{
+ Word Abp,B,B1,B2,I,I1,Ip,Ipp,L,Ls,a1,a2,ahs1,as1,b0,b1,bhs1,
+ bs1,c,ch,dp,k,m,m1,n,r,r0,r1,s,s1,s2,sbp,sbp1,sp,ss1,t0,
+ t1,tp0,tp1,ts1,u,u1,u2,up,up1,us1,v1,vp1,vs1,ws0,ws1,z,
+ z1,z2;
+ /* hide B1,B2,I1,Ls,a1,a2,b0,b1,k,m,m1,n,r,r0,r1,s,s1,s2,sbp,
+ sbp1,sp,ss1,t0,t1,tp0,tp1,ts1,ws0,ws1; */
+
+Step1: /* Ap without roots. */
+ dp = RNNEG(d);
+ if (Lp == NIL)
+ {
+ I = LIST2(dp,d);
+ L = LIST2(I,1);
+ goto Return;
+ }
+
+Step2: /* Initialize. */
+ Abp = AFUPGS(M,Ap);
+ AFUPGC(M,A,Abp,&B,&B1,&B2);
+ n = PDEG(A);
+ k = PDEG(B);
+ Ls = Lp;
+ if (EVEN(n))
+ t0 = 1;
+ else
+ t0 = -1;
+ tp0 = -t0;
+ b0 = dp;
+ I1 = FIRST(Ls);
+ a1 = FIRST(I1);
+ u1 = AFPEMV(1,M,A,AFFRN(a1));
+ s1 = AFSIGN(M,J,u1);
+ L = NIL;
+ if (t0 * s1 > 0)
+ r0 = 0;
+ else
+ {
+ r0 = 1;
+ I = LIST2(b0,a1);
+ L = COMP2(1,I,L);
+ }
+ if (EVEN(k))
+ ws0 = 1;
+ else
+ ws0 = -1;
+ m = PDEG(Abp);
+ if (EVEN(m))
+ sbp1 = 1;
+ else
+ sbp1 = -1;
+
+Step3: /* Root of A in (bi,A_{i+1}). */
+ ADV2(Ls,&I1,&m1,&Ls);
+ b1 = SECOND(I1);
+ v1 = AFPEMV(1,M,A,AFFRN(b1));
+ t1 = AFSIGN(M,J,v1);
+ if (EVEN(m1))
+ tp1 = tp0;
+ else
+ tp1 = -tp0;
+ ts1 = t1;
+ if (t1 == 0)
+ ts1 = tp1;
+ if (Ls == NIL)
+ a2 = d;
+ else
+ a2 = FIRST(FIRST(Ls));
+ if (RNCOMP(b1,a2) < 0)
+ {
+ u2 = AFPEMV(1,M,A,AFFRN(a2));
+ s2 = AFSIGN(M,J,u2);
+ if (ts1 * s2 <= 0)
+ r1 = 1;
+ else
+ r1 = 0;
+ }
+ else
+ {
+ u2 = v1;
+ s2 = t1;
+ r1 = 0;
+ }
+
+Step4: /* alpha_i a root of A. */
+ if (k == 0)
+ ws1 = ws0;
+ else
+ {
+ ws1 = AFUPSR(M,J,B,b1);
+ if (ws1 == 0)
+ ws1 = -ws0;
+ }
+ if (ws0 * ws1 < 0)
+ {
+ L = COMP2(m1 + 1,I1,L);
+ goto Step9;
+ }
+
+Step5: /* Other roots of A in (ai,bi). */
+ r = r0 + r1;
+ if (r == 2)
+ goto Step9;
+ if (s1 != 0)
+ ss1 = s1;
+ else
+ ss1 = tp0;
+ if (r == 1 || EVEN(m1))
+ {
+ if (ss1 * t1 <= 0)
+ goto Step6;
+ else
+ goto Step9;
+ }
+ if (t1 == 0)
+ {
+ if (ss1 * tp0 > 0)
+ goto Step6;
+ as1 = a1;
+ bs1 = b1;
+ goto Step8;
+ }
+ if (ss1 * t1 < 0)
+ goto Step6;
+ if (ss1 * tp0 > 0)
+ goto Step9;
+ else
+ goto Step7;
+
+Step6: /* One root in (ai,bi). */
+ I = AFPRRI(M,J,A,Abp,I1,ss1,sbp1);
+ L = COMP2(1,I,L);
+ goto Step9;
+
+Step7: /* Zero or two roots of A in (ai,bi). */
+ as1 = a1;
+ bs1 = b1;
+ us1 = u1;
+ vs1 = v1;
+ ahs1 = AFFRN(as1);
+ bhs1 = AFFRN(bs1);
+ up1 = AFPEMV(1,M,Ap,ahs1);
+ vp1 = AFPEMV(1,M,Ap,bhs1);
+ if (vp1 == 0)
+ goto Step9;
+ do
+ {
+ z = AFQ(M,us1,up1);
+ z1 = AFDIF(ahs1,z);
+ z = AFQ(M,vs1,vp1);
+ z2 = AFDIF(bhs1,z);
+ z = AFCOMP(M,J,z1,z2);
+ if (z >= 0)
+ goto Step9;
+ c = RIB(as1,bs1);
+ ch = AFFRN(c);
+ u = AFPEMV(1,M,A,ch);
+ s = AFSIGN(M,J,u);
+ up = AFPEMV(1,M,Ap,ch);
+ sp = AFSIGN(M,J,up);
+ if (ss1 * s > 0 && sp == 0)
+ goto Step9;
+ if (s == 0 && tp0 * sp < 0)
+ {
+ bs1 = c;
+ goto Step8;
+ }
+ if (ss1 * s <= 0)
+ {
+ Ip = LIST2(as1,c);
+ Ipp = LIST2(c,bs1);
+ if (tp0 * sp <= 0)
+ {
+ I = AFPRRI(M,J,A,Abp,Ip,ss1,sbp1);
+ L = COMP4(1,Ipp,1,I,L);
+ }
+ else
+ {
+ I = AFPRRI(M,J,A,Abp,Ipp,-ss1,sbp1);
+ L = COMP4(1,I,1,Ip,L);
+ }
+ goto Step9;
+ }
+ if (tp0 * sp > 0)
+ {
+ as1 = c;
+ us1 = u;
+ up1 = up;
+ }
+ else
+ {
+ bs1 = c;
+ vs1 = u;
+ vp1 = up;
+ }
+ }
+ while (1);
+
+Step8: /* Roots at bsi and in (asi,bsi) */
+ do
+ {
+ c = RIB(as1,bs1);
+ s = AFUPSR(M,J,A,c);
+ if (ss1 * s <= 0)
+ {
+ sbp = AFUPSR(M,J,Abp,c);
+ Ip = LIST2(as1,c);
+ Ipp = LIST2(c,bs1);
+ if (sbp1 * sbp <= 0)
+ {
+ I = AFPRRI(M,J,A,Abp,Ip,ss1,sbp1);
+ L = COMP4(1,Ipp,1,I,L);
+ }
+ else
+ {
+ I = AFPRRI(M,J,A,Abp,Ipp,-ss1,sbp1);
+ L = COMP4(1,I,1,Ip,L);
+ }
+ goto Step9;
+ }
+ as1 = c;
+ }
+ while (1);
+
+Step9: /* Update. */
+ if (r1 == 1)
+ {
+ I = LIST2(b1,a2);
+ L = COMP2(1,I,L);
+ }
+ a1 = a2;
+ r0 = r1;
+ tp0 = tp1;
+ s1 = s2;
+ ws0 = ws1;
+ sbp1 = -sbp1;
+ u1 = u2;
+ if (Ls != NIL)
+ goto Step3;
+
+Step10: /* Finish. */
+ L = INV(L);
+
+Return: /* Prepare for return. */
+ return(L);
+}
+@
+
+
+1.1
+log
+@Initial revision
+@
+text
+@d1 1
+@
131 RCS/AFPRLS.c,v
@@ -0,0 +1,131 @@
+head 1.2;
+access;
+symbols;
+locks
+ saclib:1.2; strict;
+comment @ * @;
+
+
+1.2
+date 95.08.09.13.56.41; author Mark; state del;
+branches;
+next 1.1;
+
+1.1
+date 95.08.09.13.44.12; author saclib; state init;
+branches;
+next ;
+
+
+desc
+@Algebraic number field polynomial real root list separation.
+@
+
+
+1.2
+log
+@Replaced by AFUPRLS.
+@
+text
+@[ removed from library ]
+/*======================================================================
+ AFPRLS(M,I,A1,A2,L1,L2; Ls1,Ls2)
+
+Algebraic number field polynomial real root list separation.
+
+Inputs
+ M : in Z[x], the minimal polynomial of an algebraic number alpha.
+ I : an acceptable isolating interval for alpha.
+ A1,A2 : in Q(alpha)[x]. A1 and A2 do not have any common roots and
+ their real roots are of odd multiplicity.
+ L1,L2 : strong isolation lists for the real roots of A1 and A2 resp.
+ L1 = (I1_1,m1_1,...,I1_r1,m1_r1),
+ L2 = (I2_1,m2_1,...,I2_r2,m2_r2).
+ I1_1 < I1_2 < ... < I1_r1 and I2_1 < I2_2 < ... < I2_r2.
+
+Outputs
+ Ls1,Ls2 : lists (Isi_1,mi_1,...,Isi*_r1,mi_r1) where Isi_j is a
+ binary rational subinterval of Ii_j containing the root
+ of Ai in Ii_j. Each Is1_j is strongly disjoint from each
+ Is2_j.
+======================================================================*/
+#include "saclib.h"
+
+void AFPRLS(M,I,A1,A2,L1,L2, Ls1_,Ls2_)
+ Word M,I,A1,A2,L1,L2, *Ls1_,*Ls2_;
+{
+ Word I1,I2,Lp1,Lp2,Ls1,Ls2,m1,m2,s;
+ /* hide Lp1,Lp2,m1,m2,s; */
+
+Step1: /* Initialize. */
+ if (L1 == NIL || L2 == NIL)
+ {
+ Ls1 = L1;
+ Ls2 = L2;
+ goto Return;
+ }
+ ADV2(L1,&I1,&m1,&Lp1);
+ Ls1 = NIL;
+ ADV2(L2,&I2,&m2,&Lp2);
+ Ls2 = NIL;
+
+Step2: /* Refine and merge. */
+ do
+ {
+ AFPRRS(M,I,A1,A2,I1,I2,&I1,&I2,&s);
+ if (s < 0)
+ {
+ Ls1 = COMP2(m1,I1,Ls1);
+ if (Lp1 != NIL)
+ ADV2(Lp1,&I1,&m1,&Lp1);
+ else
+ I1 = 0;
+ }
+ else
+ {
+ Ls2 = COMP2(m2,I2,Ls2);
+ if (Lp2 != NIL)
+ ADV2(Lp2,&I2,&m2,&Lp2);
+ else
+ I2 = 0;
+ }
+ }
+ while (!(I1 == 0 || I2 == 0));
+
+Step3: /* Finish. */
+ if (I1 == 0)
+ {
+ Ls2 = COMP2(m2,I2,Ls2);
+ while (Lp2 != NIL)
+ {
+ ADV2(Lp2,&I2,&m2,&Lp2);
+ Ls2 = COMP2(m2,I2,Ls2);
+ }
+ }
+ else
+ {
+ Ls1 = COMP2(m1,I1,Ls1);
+ while (Lp1 != NIL)
+ {
+ ADV2(Lp1,&I1,&m1,&Lp1);
+ Ls1 = COMP2(m1,I1,Ls1);
+ }
+ }
+ Ls1 = INV(Ls1);
+ Ls2 = INV(Ls2);
+
+Return: /* Prepare for return. */
+ *Ls1_ = Ls1;
+ *Ls2_ = Ls2;
+ return;
+}
+@
+
+
+1.1
+log
+@Initial revision
+@
+text
+@d1 1
+@
94 RCS/AFPRRI.c,v
@@ -0,0 +1,94 @@
+head 1.2;
+access;
+symbols;
+locks
+ saclib:1.2; strict;
+comment @ * @;
+
+
+1.2
+date 95.08.09.13.57.14; author Mark; state del;
+branches;
+next 1.1;
+
+1.1
+date 95.08.09.13.44.37; author saclib; state init;
+branches;
+next ;
+
+
+desc
+@Algebraic number field polynomial relative real root isolation.
+@
+
+
+1.2
+log
+@Replaced by AFUPRRI.
+@
+text
+@[ removed from library ]
+/*======================================================================
+ Js <- AFPRRI(M,I,A,B,J,s1,t1)
+
+Algebraic number field polynomial relative real root isolation.
+
+Inputs
+ M : in Z[x], the minimal polynomial of an algebraic number alpha.
+ I : an acceptable isolating interval for alpha.
+ A,B : in Q(alpha)[x].
+ J : a left open, right closed interval (a1,a2) where a1 and a2
+ are binary rational numbers with a1 < a2.
+ A and B have unique roots, r1 and r2 respectively, in J,
+ each of odd multiplicity and with r1 not equal r2.
+ s1,t1 : BETA-digits, s1 = sign(A(a1+)) and t1 = sign(B(a1+)).
+
+Outputs
+ Js : a left-open, right closed interval. Js = (as1,as2) is a
+ subinterval of J with as1 and as2 binary rational numbers
+ and as1 < as2, such that Js contains r1 but not r2.
+======================================================================*/
+#include "saclib.h"
+
+Word AFPRRI(M,I,A,B,J,s1,t1)
+ Word M,I,A,B,J,s1,t1;
+{
+ Word Js,a,as1,as2,s,t,u,v;
+ /* hide s,t,u,v; */
+
+Step1: /* Initialize. */
+ FIRST2(J,&as1,&as2);
+
+Step2: /* Bisect. */
+ a = RIB(as1,as2);
+ s = AFUPSR(M,I,A,a);
+ if (s == 0)
+ s = -s1;
+ t = AFUPSR(M,I,B,a);
+ if (t == 0)
+ t = -t1;
+ u = s1 * s;
+ v = t1 * t;
+ if (u > 0)
+ as1 = a;
+ else
+ as2 = a;
+ if (u == v)
+ goto Step2;
+
+Step3: /* Construct Js. */
+ Js = LIST2(as1,as2);
+
+Return: /* Prepare for return. */
+ return(Js);
+}
+@
+
+
+1.1
+log
+@Initial revision
+@
+text
+@d1 1
+@
155 RCS/AFPRRS.c,v
@@ -0,0 +1,155 @@
+head 1.2;
+access;
+symbols;
+locks
+ saclib:1.2; strict;
+comment @ * @;
+
+
+1.2
+date 95.08.09.13.53.56; author Mark; state del;
+branches;
+next 1.1;
+
+1.1
+date 95.08.09.13.42.15; author saclib; state init;
+branches;
+next ;
+
+
+desc
+@Algebraic number field polynomial real root separation.
+@
+
+
+1.2
+log
+@Replaced by AFUPRRS.
+@
+text
+@[ removed from library ]
+/*======================================================================
+ AFPRRS(M,I,A1,A2,I1,I2; Is1,Is2,s)
+
+Algebraic number field polynomial real root separation.
+
+Inputs
+ M : in Z[x], the minimal polynomial of an algebraic number alpha.
+ I : an acceptable isolating interval for alpha.
+ A1,A2 : in Q(alpha)[x], deg(A1) > 0, deg(A2) > 0.
+ I1,I2 : Binary rational intervals, each of which is either left-open
+ and right-closed, or a one-point interval. I1 contains a
+ unique root alpha1 of A1 of odd multiplicity, and I2
+ contains a unique root alpha2 not equal alpha1 of A2 of
+ odd multiplicity.
+
+Outputs
+ Is1,Is2 : Binary rational intervals. Is1 is a subinterval of I1
+ containing alpha1 and Is2 is a subinterval of I2 containing
+ alpha2. Is1 and Is2 are strongly disjoint.
+ If I1 is left-open and right-closed then so is Is1, and
+ similarly for I2 and Is2.
+ s : a BETA-digit, s = -1 if Is1 < Is2, and s = 1 if Is1 > Is2.
+======================================================================*/
+#include "saclib.h"
+
+void AFPRRS(M,I,A1,A2,I1,I2, Is1_,Is2_,s_)
+ Word M,I,A1,A2,I1,I2, *Is1_,*Is2_,*s_;
+{
+ Word Is1,Is2,a1,a2,b1,b2,c,d1,d2,s,s1,s2,t,u,v;
+ /* hide s,s1,s2,t,u,v; */
+
+Step1: /* I1 and I2 disjoint. */
+ FIRST2(I1,&a1,&b1);
+ FIRST2(I2,&a2,&b2);
+ if (RNCOMP(b1,a2) < 0)
+ {
+ Is1 = I1;
+ Is2 = I2;
+ s = -1;
+ goto Return;
+ }
+ if (RNCOMP(b2,a1) < 0)
+ {
+ Is1 = I1;
+ Is2 = I2;
+ s = 1;
+ goto Return;
+ }
+
+Step2: /* Initialize. */
+ d1 = RNDIF(b1,a1);
+ d2 = RNDIF(b2,a2);
+ s1 = 2;
+ s2 = 2;
+
+Step3: /* Bisect I1. */
+ t = RNCOMP(d1,d2);
+ if (t >= 0)
+ {
+ if (s1 > 1)
+ s1 = AFUPSR(M,I,A1,b1);
+ c = RIB(a1,b1);
+ u = AFUPSR(M,I,A1,c);
+ if (s1 == 0 || s1 * u < 0)
+ {
+ a1 = c;
+ v = 1;
+ }
+ else
+ {
+ b1 = c;
+ s1 = u;
+ v = -1;
+ }
+ d1 = RNDIF(b1,a1);
+ }
+
+Step4: /* Bisect I2. */
+ if (t < 0)
+ {
+ if (s2 > 1)
+ s2 = AFUPSR(M,I,A2,b2);
+ c = RIB(a2,b2);
+ u = AFUPSR(M,I,A2,c);
+ if (s2 == 0 || s2 * u < 0)
+ {
+ a2 = c;
+ v = -1;
+ }
+ else
+ {
+ b2 = c;
+ s2 = u;
+ v = 1;
+ }
+ d2 = RNDIF(b2,a2);
+ }
+
+Step5: /* I1 and I2 disjoint. */
+ if (v < 0 && RNCOMP(b1,a2) < 0)
+ s = -1;
+ else
+ if (v > 0 && RNCOMP(b2,a1) < 0)
+ s = 1;
+ else
+ goto Step3;
+ Is1 = LIST2(a1,b1);
+ Is2 = LIST2(a2,b2);
+
+Return: /* Prepare for return. */
+ *Is1_ = Is1;
+ *Is2_ = Is2;
+ *s_ = s;
+ return;
+}
+@
+
+
+1.1
+log
+@Initial revision
+@
+text
+@d1 1
+@
80 RCS/AFPSIP.c,v
@@ -0,0 +1,80 @@
+head 1.1;
+access;
+symbols;
+locks
+ saclib:1.1; strict;
+comment @ * @;
+
+
+1.1
+date 2000.01.01.20.20.40; author George; state init;
+branches;
+next ;
+
+
+desc
+@Algebraic field polynomial to software interval polynomial.
+@
+
+
+1.1
+log
+@Initial revision
+@
+text
+@/*======================================================================
+ AFPSIP(I,A,p,B)
+
+Algebraic field polynomial to software interval polynomial.
+
+Inputs
+ I : a standard logarithmic binary rational interval containing
+ an algebraic number \alpha.
+ A : a univariate polynomial of positive degree over the
+ algebraic number field Q(\alpha). Each nonzero coefficient
+ of A is represented by a pair (r,a) where r is a
+ rational number and a is an integral polynomial of
+ degree less than the degree of \alpha. Let n be the
+ degree of A.
+ p : a positive BETA-digit.
+ B : an array of size at least 2(n + 1)(p + 3) + 1.
+
+Effect
+ A software interval polynomial of precision p containing A(x)
+ is placed in B.
+======================================================================*/
+#include "floats.h"
+
+void AFPSIP(I,A,p,B)
+ BDigit *B,p;
+ Word I,A;
+{
+ BDigit *Ip,i,m,q1,q2,s;
+ Word a,Ap;
+
+Step1: /* Compute parameters. */
+ m = PDEG(A);
+ q1 = p + 3;
+ q2 = q1 + q1;
+ s = (m + 1) * q2 + 1;
+ B[0] = m;
+
+Step2: /* Convert isolating interval to floating point. */
+ Ip = GETARRAY(q2);
+ LBRISI(I,p,Ip);
+
+Step3: /* Convert and evaluate coefficients. */
+ Ap = A;
+ for (i = m; i >= 0; i--) {
+ s = s - q2;
+ if (Ap == NIL || PDEG(Ap) < i)
+ ISI(0,p,B + s);
+ else {
+ Ap = RED(Ap);
+ ADV(Ap,&a,&Ap);
+ ANSI(Ip,a,B + s); } }
+
+Return: /* Return. */
+ return;
+}
+@
107 RCS/AFPSIR.c,v
@@ -0,0 +1,107 @@
+head 1.2;
+access;
+symbols;
+locks
+ saclib:1.2; strict;
+comment @ * @;
+
+
+1.2
+date 2001.12.17.18.12.55; author George; state embe;
+branches;
+next 1.1;
+
+1.1
+date 2000.07.10.18.39.52; author George; state new;
+branches;
+next ;
+
+
+desc
+@Algebraic field polynomial software interval refinement.
+@
+
+
+1.2
+log
+@Fixed the header.
+@
+text
+@/*======================================================================
+ K <- AFPSIR(p,I,B,J)
+
+Algebraic field polynomial software interval refinement.
+
+Inputs
+ p : a nonnegative BETA digit (the precision).
+ I : a standard logarithmic binary rational that is a refined
+ isolating interval for an algebraic number alpha.
+ I must be exactly convertible to a precision p
+ software interval.
+ B : a univariate polynomial of positive degree over Q(alpha).
+ Each coefficient of B is represented as a rational
+ number and a univariate polynomial over Z[alpha].
+ J : a standard logarithmic binary rational isolating interval
+ for a root beta of B(x).
+
+Output
+ K : a standard logarithmic binary rational isolating interval
+ for beta, as small as could be achieved using precision p
+ software interval arithmetic, where p is the precision of I.
+======================================================================*/
+#include "saclib.h"
+
+Word AFPSIR(p,I,B,J)
+ BDigit p;
+ Word I,B,J;
+{
+ BDigit *Bp,n,q,S,t,*T1,*T2;
+ Word b2,J2,K;
+
+Step1: /* Convert B to a software interval polynomial. */
+ q = p + p + 6;
+ n = PDEG(B);
+ S = (n + 1) * q + 1;
+ Bp = GETARRAY(S);
+ AFPSIP(I,B,p,Bp);
+
+Step2: /* Compute the trend of beta. */
+ b2 = SECOND(J);
+ J2 = LIST2(b2,b2);
+ T1 = GETARRAY(q);
+ LBRISI(J2,p,T1);
+ T2 = GETARRAY(q);
+ SIPEVAL(Bp,T1,T2);
+ FREEARRAY(T1);
+ t = SISIGN(T2);
+ FREEARRAY(T2);
+ if (t == NIL) {
+ K = J;
+ goto Step4; }
+
+Step3: /* Refine J. */
+ K = SIPIR(Bp,J,t,-p * ZETA);
+
+Step4: /* Free array Bp. */
+ FREEARRAY(Bp);
+
+Return: /* Return K. */
+ return(K);
+}
+@
+
+
+1.1
+log
+@Initial revision
+@
+text
+@d2 1
+a2 1
+ K = AFPSIR(p,I,B,J)
+d10 1
+a10 1
+ I must be exactly convertivle to a precision p
+a23 1
+#include "dcell.h"
+@
187 RCS/AFPSUM.c,v
@@ -0,0 +1,187 @@
+head 1.2;
+access;
+symbols;
+locks
+ saclib:1.2; strict;
+comment @ * @;
+
+
+1.2
+date 94.05.11.14.27.43; author Mark; state spec;
+branches;
+next 1.1;
+
+1.1
+date 94.05.11.14.25.14; author saclib; state init;
+branches;
+next ;
+
+
+desc
+@Algebraic number field polynomial sum.
+@
+
+
+1.2
+log
+@Input r = 0 is OK.
+@
+text
+@/*===========================================================================
+ C <- AFPSUM(r,A,B)
+
+Algebraic number field polynomial sum.
+
+Inputs
+ r : a non-negative BETA-digit.
+ A,B : in Q(alpha)[x_1,...,x_r].
+
+Outputs
+ C : in Q(alpha)[x_1,...,x_r], C = A + B.
+===========================================================================*/
+#include "saclib.h"
+
+Word AFPSUM(r,A,B)
+ Word r,A,B;
+{
+ Word Ap,Bp,C,Cp,a,b,c,e,f,rp;
+
+Step1: /* A or B zero. */
+ if (A == 0) {
+ C = B;
+ goto Return; }
+ if (B == 0) {
+ C = A;
+ goto Return; }
+
+Step2: /* r = 0. */
+ if (r == 0) {
+ C = AFSUM(A,B);
+ goto Return; }
+
+Step3: /* Match coefficients. */
+ Ap = A;
+ Bp = B;
+ Cp = NIL;
+ rp = r - 1;
+ do {
+ e = FIRST(Ap);
+ f = FIRST(Bp);
+ if (e > f) {
+ ADV2(Ap,&e,&a,&Ap);
+ Cp = COMP2(a,e,Cp); }
+ else
+ if (e < f) {
+ ADV2(Bp,&f,&b,&Bp);
+ Cp = COMP2(b,f,Cp); }
+ else {
+ ADV2(Ap,&e,&a,&Ap);
+ ADV2(Bp,&f,&b,&Bp);
+ if (rp == 0)
+ c = AFSUM(a,b);
+ else
+ c = AFPSUM(rp,a,b);
+ if (c != 0)
+ Cp = COMP2(c,e,Cp); } }
+ while (Ap != NIL && Bp != NIL);
+
+Step4: /* Finish. */
+ if (Ap == NIL)
+ Ap = Bp;
+ if (Cp == NIL)
+ C = Ap;
+ else {
+ C = INV(Cp);
+ SRED(Cp,Ap); }
+ if (C == NIL)
+ C = 0;
+
+Return: /* Prepare for return. */
+ return(C);
+}
+@
+
+
+1.1
+log
+@Initial revision
+@
+text
+@d1 2
+a2 2
+/*======================================================================
+ C <- AFPSUM(r,A,B)
+d7 2
+a8 2
+ r : a BETA-digit, r >= 1.