{VERSION 6 0 "SUN SPARC SOLARIS" "6.0" } {USTYLETAB {CSTYLE "Maple Input" -1 0 "Courier" 0 1 255 0 0 1 0 1 0 0 1 0 0 0 0 1 }{CSTYLE "" -1 256 "" 1 24 0 0 0 0 0 0 0 0 0 0 0 0 0 1 } {CSTYLE "" -1 257 "" 0 1 0 0 0 0 0 1 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 258 "" 1 24 0 0 0 0 0 0 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 259 "" 0 1 0 0 0 0 0 1 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 260 "" 0 1 0 0 0 0 0 1 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 261 "" 0 1 0 0 0 0 0 1 0 0 0 0 0 0 0 1 } {CSTYLE "" -1 262 "" 0 1 0 0 0 0 0 1 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 263 "" 0 1 0 0 0 0 0 1 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 264 "" 0 1 0 0 0 0 0 1 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 265 "" 0 1 0 0 0 0 0 1 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 266 "" 0 1 0 0 0 0 0 1 0 0 0 0 0 0 0 1 } {CSTYLE "" -1 267 "" 0 1 0 0 0 0 0 1 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 268 "" 0 1 0 0 0 0 0 1 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 269 "" 0 1 0 0 0 0 0 1 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 270 "" 0 1 0 0 0 0 0 1 0 0 0 0 0 0 0 1 }{CSTYLE "" -1 271 "" 0 1 0 0 0 0 0 1 0 0 0 0 0 0 0 1 } {PSTYLE "Normal" -1 0 1 {CSTYLE "" -1 -1 "" 0 1 0 0 0 0 0 0 0 0 0 0 0 0 0 1 }0 0 0 -1 -1 -1 0 0 0 0 0 0 -1 0 }{PSTYLE "" 0 256 1 {CSTYLE "" -1 -1 "" 0 1 0 0 0 0 0 1 0 0 0 0 0 0 0 1 }0 0 0 -1 -1 -1 0 0 0 0 0 0 -1 0 }{PSTYLE "" 0 257 1 {CSTYLE "" -1 -1 "" 0 1 0 0 0 0 0 1 0 0 0 0 0 0 0 1 }0 0 0 -1 -1 -1 0 0 0 0 0 0 -1 0 }{PSTYLE "" 0 258 1 {CSTYLE " " -1 -1 "" 0 1 0 0 0 0 0 1 0 0 0 0 0 0 0 1 }0 0 0 -1 -1 -1 0 0 0 0 0 0 -1 0 }} {SECT 0 {EXCHG {PARA 0 "" 0 "" {TEXT 256 30 "Greatest Common Divisor, \+ gcd." }}{PARA 0 "" 0 "" {TEXT -1 0 "" }}{PARA 0 "" 0 "" {TEXT -1 55 " In Maple: gcd for polynomials, and igcd for integers.\n " }}{PARA 0 " " 0 "" {TEXT 264 9 "Notation:" }{TEXT -1 54 " Z = set of all integers = \{...,-3,-2,-1,0,1,2,3,...\}" }}{PARA 0 "" 0 "" {TEXT -1 0 "" }} {PARA 0 "" 0 "" {TEXT -1 4 "The " }{TEXT 263 6 "ideal " }{TEXT -1 106 "(a1, a2, a3, .., ak) is defined as:\n Z*a1 + Z*a2 + ... + Z*ak = \{n1 *a1+n2*a2+...+nk*ak where n1..nk in Z\}" }}{PARA 0 "" 0 "" {TEXT -1 0 "" }}{PARA 0 "" 0 "" {TEXT -1 87 "Note: (a1,a2,a3,...,ak) = (a2,a3, ..ak) if and only if a1 is an element of (a2,..,ak)." }}{PARA 0 "" 0 "" {TEXT -1 111 "More specifically: (a1, a2) = (a1) if and only if a 2 is in (a1), if and only if a2 = n*a1 for some integer n." }}{PARA 0 "" 0 "" {TEXT -1 47 "Furthermore:\n(a1-a2, a2) = (a1,a2) (rule 1) " }}{PARA 0 "" 0 "" {TEXT -1 3 "and" }}{PARA 0 "" 0 "" {TEXT -1 38 "(a 2, a1) = (a1, a2). (rule 2)" }}{PARA 0 "" 0 "" {TEXT -1 3 "an d" }}{PARA 0 "" 0 "" {TEXT -1 42 "(0, a1) = (a1). (r ule 3)" }}{PARA 0 "" 0 "" {TEXT -1 0 "" }}{PARA 0 "" 0 "" {TEXT -1 38 "We will show that rule (1) is true by " }{TEXT 257 11 "an example:" } }{PARA 0 "" 0 "" {TEXT -1 37 " (14, 8) = \{n*14 + m*8 | n,m in Z\} " }}{PARA 0 "" 0 "" {TEXT -1 187 " = \{n*(14-8) + (m+ n)*8 | n,m in Z\}\nand now if we take all m's in Z, then that's the sa me as taking all m+n's in Z, so the latter equals:\n \{n*(14-8) + m*8 | n,m in Z\}" }}{PARA 0 "" 0 "" {TEXT -1 26 " = \{n*6 + m*8 | n, m in Z\}" }}{PARA 0 "" 0 "" {TEXT -1 42 " = (6,8)\n\nUsing rules (1),( 2),(3) we find:" }}{PARA 0 "" 0 "" {TEXT -1 120 "(14,8) = (14-8,8) = ( 6,8) = (8,6) = (8-6,6) = (2,6) = (6,2) = (6-2,2) = (4,2) = (4-2,2) = ( 2,2) = (2-2,2) = (0,2) = (2)." }}{PARA 0 "" 0 "" {TEXT -1 0 "" }} {PARA 256 "" 0 "" {TEXT -1 130 "Definition: if a1..ak are integers, th en the gcd of a1..ak is defined as the non-negative integer d for whic h (a1,a2,..,ak) = (d)." }}{PARA 0 "" 0 "" {TEXT -1 0 "" }}{PARA 0 "" 0 "" {TEXT -1 261 "Remark: Since a1 is an element of (a1..ak) and thus an element of (d) = \{n*d | n in Z\} we see that there exists an inte ger n such that a1 = n*d. So a1 is a multiple of d, and d is a divisor of a1. Likewise, d is also a divisor of each a.i for i=1..k, so d is \+ a " }{TEXT 259 14 "common divisor" }{TEXT -1 11 " of a1..ak." }}{PARA 0 "" 0 "" {TEXT -1 104 "Conversely d is also an element of (d) = (a1,. .,ak) = \{n1*a1 + .. + nk*ak where n1..nk in Z\}. Therefore:" }}{PARA 0 "" 0 "" {TEXT 260 38 "there exist integers n1..nk such that:" }} {PARA 0 "" 0 "" {TEXT -1 30 " d = n1*a1 + ... + nk*ak." }}{PARA 0 "" 0 "" {TEXT -1 0 "" }}{PARA 0 "" 0 "" {TEXT -1 161 "In our example we can find for example 2 = (-1)*14 + 2*8 but there are also other \+ solutions, for example: 2 = 3*14 - 5*8. So we see that n1..nk are not \+ unique." }}{PARA 0 "" 0 "" {TEXT -1 130 "Note that if g divides a1..ak (i.e. each a1..ak is a multiple of g) then g also divides n1*a1+...+n k*ak and hence g divides d. So:" }}{PARA 0 "" 0 "" {TEXT 261 52 "Any i nteger g that divides a1..ak divides d as well." }{TEXT -1 92 " Theref ore, we can say that d is not only a common divisor of a1,..,ak, but w e can say that " }{TEXT 262 32 "d is the greatest common divisor" } {TEXT -1 11 " of a1..ak." }}{PARA 0 "" 0 "" {TEXT -1 0 "" }}{PARA 0 " " 0 "" {TEXT -1 84 "Lets calculate the gcd of two nonnegative integers a,b using rules (1), (2) and (3)." }}{PARA 0 "" 0 "" {TEXT -1 277 "Re member we are searching for an integer d such that (a,b) = (d). We do \+ this by decreasing a,b with rule (1) while keeping the ideal (a,b) th e same. Since a,b are nonnegative and they keep decreasing, at some po int one of them must be 0, then the other is the gcd by rule (3)." }} {PARA 0 "" 0 "" {TEXT -1 0 "" }}{PARA 0 "" 0 "" {TEXT -1 200 "Note: \+ if you do proc(a::nonnegint, b::nonnegint) instead of proc(a,b) then t he procedure will only accept nonnegative integers as input, and will \+ give an error message on inputs of a different type." }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 234 "GCD := proc(a::nonnegint, b::nonnegint) \n if a=0 then\n b\n elif b=0 then\n a\n elif a>b th en\n # rule (1):\n procname(a-b,b)\n else\n # rules ( 1)+(2) give: (a,b-a) = (a,b)\n procname(a,b-a)\n fi\nend;\n" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 13 "GCD( 12, 18);" }}}{EXCHG {PARA 0 "" 0 "" {TEXT -1 72 "Now the same algorithm, but a little bit \+ faster, based on the fact that " }{TEXT 265 48 "the ideal (a,b) equals the ideal (b, irem(a,b))." }}{PARA 0 "" 0 "" {TEXT -1 0 "" }}{PARA 0 "" 0 "" {TEXT -1 24 "irem = integer remainder" }}{PARA 0 "" 0 "" {TEXT -1 24 "iquo = integer quotient." }}{PARA 0 "" 0 "" {TEXT -1 0 " " }}{PARA 0 "" 0 "" {TEXT -1 102 "If a,b are integers and b>0 then we \+ can do a long division (learned in elementary school) and we find:" }} {PARA 0 "" 0 "" {TEXT -1 13 " a = q*b + r" }}{PARA 0 "" 0 "" {TEXT -1 50 "where q is the quotient, in Maple: q = iquo(a,b)," }}{PARA 0 " " 0 "" {TEXT -1 48 "and r is the remainder, in Maple: r = irem(a,b)." }}{PARA 0 "" 0 "" {TEXT -1 43 "Furthermore we have that r in \{0,1,.., b-1\}." }}{PARA 0 "" 0 "" {TEXT -1 0 "" }}{PARA 0 "" 0 "" {TEXT -1 32 "For example: if a=15, b=6 then:" }}{PARA 0 "" 0 "" {TEXT -1 58 "15 = 2*6 + 3 so the quotient is 2 and the remainder is 3." }}{PARA 0 "" 0 "" {TEXT -1 74 "The remainder is 0 if and only if b divides a, so if a is a multiple of b." }}{PARA 0 "" 0 "" {TEXT -1 0 "" }}{PARA 0 "" 0 "" {TEXT -1 26 "Now if q is an integer and" }}{PARA 0 "" 0 "" {TEXT -1 12 " a = q*b + r" }}{PARA 0 "" 0 "" {TEXT -1 5 "then:" }}{PARA 0 " " 0 "" {TEXT -1 45 "(a,b) = \{n*a + m*b where n,m = all integers\}" } }{PARA 0 "" 0 "" {TEXT -1 54 " = \{n*r + (m+n*q)*b where n,m \+ = all integers\}" }}{PARA 0 "" 0 "" {TEXT -1 47 " = \{n*r + m*b where n,m = all integers\}" }}{PARA 0 "" 0 "" {TEXT -1 23 " = (r,b) = (b,r)" }}{PARA 0 "" 0 "" {TEXT -1 25 " = (b, irem(a,b) )." }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 434 "GCD := proc(a::nonne gint, b::nonnegint)\n if b=0 then\n a\n else\n # Now b<> 0 so I can replace a by irem(a,b)\n # without changing the ideal \+ (a,b) = Z*a+Z*b.\n # I'll also interchange a and b, this way from \n # this point on, b will always be smaller than a,\n # whi ch causes that the irem in the next step will\n # always make a s maller, so the algorithm will end.\n procname(b,irem(a,b))\n fi \nend;" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 17 "GCD(100!, 12345); " }}}{EXCHG {PARA 0 "" 0 "" {TEXT -1 172 "Now lets do the same for pol ynomials. If f,g in Q[x] then we will compute the gcd d. Recall that Q [x] stands for all polynomials in x with rational numbers as coefficie nts." }}{PARA 0 "" 0 "" {TEXT -1 0 "" }}{PARA 0 "" 0 "" {TEXT -1 78 "A gain, the gcd d is defined as a polynomial d in Q[x] such that: (f,g) = (d)." }}{PARA 0 "" 0 "" {TEXT -1 60 "where the ideal (f1..fn) is de fined as: Q[x]*f1+...+Q[x]*fn." }}{PARA 0 "" 0 "" {TEXT -1 0 "" }} {PARA 0 "" 0 "" {TEXT -1 416 "Note that if c is a non-zero rational nu mber, then the ideal (d) equals the ideal (c*d). Therefore, if d is a correct answer of the following procedure, then c*d is also a correct answer for any non-zero rational number c. In a later procedure belo w we will divide the gcd d by a suitable rational number c (the conten t of the polynomial) in order to obtain a gcd that is smaller (i.e. oc cupies less screen space)." }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 133 "GCD:=proc(f,g,x)\n if g=0 then\n f\n else\n # Note: args = all arguments = f,g,x\n procname(g, rem(args), x)\n fi \nend;" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 0 "" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 29 "f:=expand( mul(x-i,i=-2..6));" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 30 "g:=expand( mul(x-i,i=4..10) );" }}} {EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 11 "GCD(f,g,x);" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 13 "content(%,x);" }}}{EXCHG {PARA 0 " " 0 "" {TEXT -1 4 "The " }{TEXT 266 8 "content " }{TEXT -1 31 "c of a \+ polynomial F in Q[x] is " }{TEXT 267 52 "the largest rational number c such that F/c in Z[x]." }{TEXT -1 300 " If you divide the gcd d by th e content c, then you get a new gcd d/c but this time all the coeffici ents are integers, and furthermore these integers are as small as poss ible (these integers have gcd = 1, so we can not divide all coefficien ts by some integer > 1 and still keep integer coefficients)." }}} {EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 7 "%% / %;" }}}{EXCHG {PARA 0 " " 0 "" {TEXT -1 266 "The polynomial F/content(F,x) is called the primi tive part of F. It is called primpart in Maple. We can use primpart to get rid of the fractions, and we can also use it to remove common int eger factors of the coefficients. In general we can use it for the fol lowing:" }}{PARA 0 "" 0 "" {TEXT -1 0 "" }}{PARA 0 "" 0 "" {TEXT -1 2 " " }{TEXT 268 181 "primpart divides the polynomial by some rational \+ number, namely by content(F,x), in such a way that the resulting polyn omial is \"as simple as possible\" in the sense explained above." }}} {EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 11 "GCD(f,g,x);" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 12 "primpart(%);" }}}{EXCHG {PARA 0 "> \+ " 0 "" {MPLTEXT 1 0 27 "primpart( 15*x+27*x^2 , x);" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 24 "f:=randpoly(x,degree=3);" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 53 "f,g := seq(expand(f*randpoly(x,degr ee=15)) , i=1..2);" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 11 "GCD(f ,g,x);" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 12 "primpart(%);" }}} {EXCHG {PARA 0 "" 0 "" {TEXT -1 314 "As you can see, the output that o ur GCD produces is not very nice, although in principle it is correct \+ because it is a generator of the Q[x]-ideal (f,g) and so we can call i t the gcd of f and g. Note that this gcd is only determined up to a no n-zero factor in Q. We could get rid of these large numbers as follows :" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 108 "GCD:=proc(f,g,x::name )\n if g=0 then\n primpart(f,x)\n else\n procname(g, rem (args), x)\n fi\nend;" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 11 " GCD(f,g,x);" }}}{EXCHG {PARA 0 "" 0 "" {TEXT -1 323 "Now it looks nice r, but it's still not very good. This algorithm will be very slow on l arger inputs because it internally computes with very large rational n umbers. Only and the end do these large rational numbers get reduced t o small integers because only then the \"primpart\" command is called \+ in the above implementation." }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 24 "f:=randpoly(x,degree=3);" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 53 "f,g := seq(expand(f*randpoly(x,degree=40)) , i=1..2);" }}} {EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 21 "start_time := time();" }}} {EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 11 "GCD(f,g,x);" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 19 "end_time := time();" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 22 "end_time - start_time;" }}}{EXCHG {PARA 0 "" 0 "" {TEXT -1 106 "That already took quite a while. For hig her degrees it will get much worse. Lets see why it takes so long:" }} }{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 136 "GCD:=proc(f,g,x::name)\n \+ lprint( length( [args] ));\n if g=0 then\n primpart(f)\n els e\n procname(g, rem(args), x)\n fi\nend;" }}}{EXCHG {PARA 0 "> \+ " 0 "" {MPLTEXT 1 0 11 "GCD(f,g,x);" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 18 "length( [f,g,x] );" }}}{EXCHG {PARA 0 "" 0 "" {TEXT -1 154 "So the sum of the number of digits of the two polynomials went up to about 70000, even though in the input its only 400 something, a nd in the output it's:" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 11 "l ength(%%);" }}}{EXCHG {PARA 0 "" 0 "" {TEXT -1 236 "So our GCD algorit hm is very slow, in fact, it's exponentially slow. Try to run some exa mples for degree 50, 100, 150, 200,... and see what happens to the num ber of digits. Lets also just look at the largest rational number that occurs." }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 207 "GCD:=proc(f,g, x::name)\n local i;\n lprint( [length( [args] ),\n max(seq(len gth(i),i=[coeffs(f,x),coeffs(g,x)]))] );\n if g=0 then\n prim part(f)\n else\n procname(g, rem(args), x)\n fi\nend;" }}} {EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 11 "GCD(f,g,x);" }}}{EXCHG {PARA 0 "" 0 "" {TEXT -1 349 "We shouldn't use \"options trace\" in th is procedure, the numbers are so large that when you'd print them on t he screen it would be a big mess. Some of those coefficients had 5000 \+ digits (in numerator and denominator combined). Some reduction during \+ the computation will help keeping coefficients smaller, and so that sh ould make the algorithm faster." }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 219 "GCD:=proc(f,g,x::name)\n local i;\n lprint( [length( [arg s] ),\n max(seq(length(i),i=[coeffs(f,x),coeffs(g,x)]))] );\n \+ if g=0 then\n primpart(f)\n else\n procname(g, primpart(re m(args),x), x)\n fi\nend;" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 11 "GCD(f,g,x);" }}}{EXCHG {PARA 0 "" 0 "" {TEXT -1 286 "Now we're loo king at close to 150 digits, which is better, but still not very good \+ when you consider that f, g and the gcd all had small coefficients. Ho wever, just adding this primpart did already speed up the code. The fo llowing example would take *much* longer without that primpart." }}} {EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 24 "f:=randpoly(x,degree=4);" }} }{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 62 "f,g := seq(expand(f*randpol y(x,degree=50,terms=10)) , i=1..2);" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 11 "GCD(f,g,x);" }}}{EXCHG {PARA 0 "" 0 "" {TEXT -1 276 " Now lets see how long it takes in Maple. Note that in maple's gcd comm and you should *NOT* give the third variable x, that leads to strange \+ errors. Maple's gcd is in fact for the ring Q[several variables] and i t will figure out by itself what the names of those variables are." }} }{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 9 "gcd(f,g);" }}}{EXCHG {PARA 0 "" 0 "" {TEXT -1 149 "So this goes pretty much instantly. Why is thi s so much faster? Optimized and compiled code or something like that? \+ If you think that, then try this:" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 32 "f:=randpoly(x,degree=250,dense):" }}}{EXCHG {PARA 0 " > " 0 "" {MPLTEXT 1 0 54 "f,g := seq(expand(f*randpoly(x,degree=250)) \+ , i=1..2):" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 0 "" }}}{EXCHG {PARA 0 "" 0 "" {TEXT -1 482 "If we applied our algorithm GCD on these f and g, which have degree 500, then we would make the computer work \+ with integers that have many million digits. Even written in machine l anguage or C and highly optimized, that would still be an extremely le nghty computation. And besides, if you have say 128 Megabyte, you can \+ not store integers of more than 128*10^6*log(256)/log(10) digits so ev en if the CPU was infinitely fast the computation would still fail, yo u'd run out of memory." }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 12 "d :=gcd(f,g):" }}}{EXCHG {PARA 0 "" 0 "" {TEXT -1 288 "As we see, the co mputation still takes less than a second. If we had computed over the \+ ring Z[x], which our algorithm GCD does, it would have been impossible to do this, to handle the many-millions of digits in such a short tim e. So we must conclude that Maple uses a different algorithm." }}} {EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 110 "GCD_mod_p:=proc(f,g,x,p)\n \+ if g=0 then\n f\n else\n procname(g, Rem(f,g,x) mod p, x, p)\n fi\nend;" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 25 "d13: =GCD_mod_p(f,g,x,13);" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 0 "" } }}{EXCHG {PARA 0 "" 0 "" {TEXT -1 575 "As we see, this algorithm GCD_m od_p is *much* faster than GCD, while the only difference is that we r educe coefficients modulo the prime number p. This makes a big differe nce because by reducing modulo p our coefficients can never grow large , and whether we use a simple multiplication algorithm or Karatsuba or even more advanced algorithms for integer arithmetic, it's obvious th at computing with just a few digits is always much faster than with ve ry large integers. We will see later in this course how gcd's can be c omputed by computing a few gcd's modulo prime numbers." }}{PARA 0 "" 0 "" {TEXT -1 0 "" }}{PARA 0 "" 0 "" {TEXT -1 0 "" }}{PARA 0 "" 0 "" {TEXT 258 24 "Computing the remainder." }}{PARA 0 "" 0 "" {TEXT -1 0 " " }}{PARA 0 "" 0 "" {TEXT -1 113 "Let's look a little closer at this r emainder procedure. How to compute the quotient and remainder of polyn omials?" }}{PARA 0 "" 0 "" {TEXT -1 911 " Basically, the idea is exact ly the same as the integer division. When you do an integer division o f positive integers a and b, then as long as you can substract a multi ple of b from a you do so. Then what you end up with is the remainder \+ of a modulo b. And, if you keep track of which multiple of b you subst racted from a (which you do in the long division algorithm that you l earned as a child), then you have the quotient too. We will see here t hat for polynomials it works exactly the same. If F,G are polynomials \+ and I want to divide F by G, then as long as the degree of F is >= the degree of G, I can substract some multiple of G from F in such a way \+ that the degree gets smaller. Note that making the degree smaller is e asy: you just have to make sure that the leading coefficient, lcoeff(F ,x) becomes 0. We have to figure out which multiple of G to subtract f rom F in order to make lcoeff(F,x) vanish." }}}{EXCHG {PARA 0 "> " 0 " " {MPLTEXT 1 0 17 "F:=3*x^5+2*x^2+1;" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 13 "G:=2*x^2+x+1;" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 59 "a := lcoeff(F,x)/lcoeff(G,x) * x^(degree(F,x)-degree(G,x));" } }}{EXCHG {PARA 0 "" 0 "" {TEXT -1 133 "Now we see that we can substrac t a multiple of G from F in such a way that the degree of G will get s maller, namely we substract a*G." }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 17 "F:=expand(F-a*G);" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 59 "a := lcoeff(F,x)/lcoeff(G,x) * x^(degree(F,x)-degree( G,x));" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 17 "F:=expand(F-a*G); " }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 59 "a := lcoeff(F,x)/lcoeff (G,x) * x^(degree(F,x)-degree(G,x));" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 17 "F:=expand(F-a*G);" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 59 "a := lcoeff(F,x)/lcoeff(G,x) * x^(degree(F,x)-degree( G,x));" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 17 "F:=expand(F-a*G); " }}}{EXCHG {PARA 0 "" 0 "" {TEXT -1 254 "Now the degree of F is small er than the degree of G. So we can no further reduce F modulo G, we ha ve found the remainder of the division of F by G. Of course, the quoti ent is just the sum of all the a's we computed. This leads to the foll owing algorithm." }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 265 "quo_re m:=proc(f,g,x)\n local F,q,a;\n F:=f;\n q:=0;\n while degree(F ,x) >= degree(g,x) do\n a := lcoeff(F,x)/lcoeff(g,x) * x^(degree(F ,x)-degree(g,x));\n q := q+a;\n F:=expand(F-a*g)\n od;\n # Now F is the remainder and q is the quotient.\n q,F\nend;" }}} {EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 14 "f:=x^12+x^4+1;" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 11 "g:=2*x^3+3;" }}}{EXCHG {PARA 0 "> \+ " 0 "" {MPLTEXT 1 0 22 "q, r :=quo_rem(f,g,x);" }}}{EXCHG {PARA 0 "> \+ " 0 "" {MPLTEXT 1 0 6 "q*g+r;" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 10 "expand(%);" }}}{EXCHG {PARA 0 "" 0 "" {TEXT -1 51 "In Maple, the syntax for quotient and remainder is:" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 21 "q := quo(f,g,x, 'r');" }}}{EXCHG {PARA 0 "" 0 "" {TEXT -1 201 "The command quo returns the quotient, which ends up in q this way. The argument 'r' is optional, and if you use that argument \+ then Maple will put the remainder in the variable r. Lets check the re sult:" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 21 "expand( q*g + r) = f;" }}}{EXCHG {PARA 0 "" 0 "" {TEXT -1 51 "For the algorithm to be co rrect we need two things:" }}{PARA 0 "" 0 "" {TEXT -1 145 "1) degree( r,x) < degree(g,x).\n2) f = q*g + r.\nBoth hold, so the algorithm is c orrect. Note that q, r are uniquely defined by conditions 1), 2)." }} }{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 15 "f:=randpoly(x);" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 6 "c:=12;" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 7 "g:=x-c;" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 12 " subs(x=c,f);" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 11 "rem(f,g,x); " }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 25 "subs(x=c, f^3+x^3*f+3 \+ );" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 29 "rem( f^3+x^3*f+3 ,g,x);" }}}{EXCHG {PARA 0 "" 0 "" {TEXT -1 273 "We see that taking t he remainder modulo x-c is the same thing as substituting x=c.\nNow ho w about taking the remainder modulo a polynomial of degree > 1?\nLet's take the polynomial g:=x^2+1\nIn Maple the complex number i is denote d by I. So I is one of the two roots of x^2+1." }}}{EXCHG {PARA 0 "> \+ " 0 "" {MPLTEXT 1 0 2 "f;" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 12 "subs(x=I,f);" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 9 "g:=x^2+1 ;" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 11 "rem(f,g,x);" }}} {EXCHG {PARA 0 "" 0 "" {TEXT -1 78 "Now another example. Lets take the number 2^(1/3). This is a root of g:=x^3-2." }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 25 "f:=randpoly(x,degree=10);" }}}{EXCHG {PARA 0 "> \+ " 0 "" {MPLTEXT 1 0 11 "c:=2^(1/3);" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 9 "g:=x^3-2;" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 12 "subs(x=c,f);" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 11 "rem(f,g ,x);" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 12 "subs(x=c,%);" }}} {EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 0 "" }}}{EXCHG {PARA 0 "" 0 "" {TEXT -1 33 "The conclusion is the following.\n" }}{PARA 257 "" 0 "" {TEXT -1 107 "If f and g are polynomials, and r:=rem(f,g,x), and if c \+ is a root of g, then: subs(x=c,f) = subs(x=c, r)." }}{PARA 0 "" 0 " " {TEXT -1 48 "This is very important, because this shows that " } {TEXT 269 41 "one can compute with an algebraic number " }{TEXT -1 39 "(i.e. with a solution of a polynomial) " }{TEXT 270 32 "by doing rema inder computations." }}{PARA 0 "" 0 "" {TEXT -1 0 "" }}{PARA 0 "" 0 " " {TEXT -1 106 "The explanation is easy,\n f = q*g+r, subs(x=c,f) = subs(x=c,q)*subs(x=c,g) + subs(x=c,r) = subs(x=c,r)." }}{PARA 0 "" 0 "" {TEXT -1 0 "" }}{PARA 0 "" 0 "" {TEXT -1 198 "So we see that comput ing in the ring Q[ 2^(1/3) ], the set of all sums of rational multiple s of powers of 2^(1/3), is completely equivalent to computing in the r ing Q[x] modulo the polynomial x^3-2." }}{PARA 258 "" 0 "" {TEXT -1 229 "Computing modulo x^3-2 simply means that whenever you see a polyn omial f with degree(f,x) >= degree(x^3-2,x) then you just replace f wi th the remainder of f modulo x^3-2, which you can calculate in Maple b y typing rem(f,x^3-2,x)." }}{PARA 0 "" 0 "" {TEXT -1 216 "In fact this proves that the ring Q[2^(1/3)] (the rational numbers extended with t he number 2^(1/3))\nis isomorphic to the ring Q[x]/(x^3-2) (polynomial s in the variable x, but we compute modulo the polynomial x^3-2)." }} {PARA 0 "" 0 "" {TEXT -1 177 "Since x^3-2 is irreducible we know from \+ algebra that these two rings are in fact fields, which means that one \+ can divide by any non-zero element. For divisions we will need the " } {TEXT 271 28 "Extended Euclidean Algorithm" }{TEXT -1 26 ", from the n ext worksheet." }}{PARA 0 "" 0 "" {TEXT -1 0 "" }}{PARA 0 "" 0 "" {TEXT -1 702 "Let m be a positive integer. Previously we have already \+ seen the ring Z/(m), the integers modulo the number m. Computing modul o m meant that if a number was negative, we add a multiple of m to mak e it non-negative, and if the number is >=m we substract multiples of \+ m to make it smaller than m. This is done by the command mod. If a is \+ an integer then \"a mod m\" will produce the remainder of a modulo m, \+ which is the same as irem(a,m), and the result is in the range 0..m-1. So every element of Z/(m) can be represented uniquely by a number in \+ this range. So the ring Z/(m) contains m elements, because in Z/(m) we will only consider integers to be different when their remainder modu lo m is different." }}{PARA 0 "" 0 "" {TEXT -1 442 "From algebra we kn ow that Z/(m) is a field (meaning that we can divide by non-zero eleme nts) if and only if m is irreducible, so m has to be a prime number p. When m is a prime number p we will denote this field by Fp. So F5 con tains five elements that we can represent by 0,1,2,3,4. The product 2* 3 for example will be 1 because 2*3 is 1+5 which reduces to 1 modulo 5 . Because 2*3 is 1 we conclude that 1/2 will be 3 and that 1/3 will be 2. " }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 10 "1/2 mod 5;" }}} {EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 0 "" }}}}{MARK "0 2 0" 55 } {VIEWOPTS 1 1 0 1 1 1803 1 1 1 1 }{PAGENUMBERS 0 1 2 33 1 1 }