A14 ==== I'm only a little embarassed for not having realized what was happening. (Perhaps I should have slept on it.) Surely I'm not alone in thinking this symbolism is highly nonintuitive. But of course, for it to be otherwise would require another protected symbol... --- Selwyn >>In: DSolve[y*D[u[x, y],x] == x*D[u[x, y],y], u[x,y], {x, y}] >> >>Out: {{u[x, y] -> C[1][(1/2)*(x^2 + y^2)]}} >> >>Square brackets are used as grouping symbols in the result!?? :^O >> >>Somebody say it isn't so. >> > It isn't so The square bracket is not delineating a factor it is enclosing the argument > to an arbitrary function named C[1]. While the function is dependent on both > x and y the dependence only occurs in the combination (x^2+y^2). > Bob Hanlon > Reply-To: ==== It isn't so. The solution is an arbitrary function of (1/2)*(x^2 + y^2)]}}. Bobby -----Original Message----- ==== Does anyone know what happened to the < in Mathematica 4.x? In some version I know I used it, but it seems to > have gone away. It allowed for real time manipulation of 3D graphics. > Ray Gittings ==== it's still there << RealTime3D` but MathGL3d may be the better solution > http://phong.informatik.uni-leipzig.de/~kuska/mathgl3dv3/ Jens >>Does anyone know what happened to the <>in Mathematica 4.x? In some version I know I used it, but it seems to >>have gone away. It allowed for real time manipulation of 3D graphics. >> >> >>Ray Gittings ==== The more I play with the example the more > depressing it gets. Start > with floating point numbers but explicitly > arbitrary-precision ones. In[1]:= > a=77617.00000000000000000000000000000; > b=33095.00000000000000000000000000000; In[3]:= > !(333.7500000000000000000000000000000 b^6 + > a^2 ((11 a^2 > b^2 - > b^6 - 121 b^4 - 2)) + > 5.500000000000000000000000000000 b^8 + > a/(2 > b)) Out[3]= > > !((-4.78339168666055402578083604864320577443814`26.6715*^32)) In[4]:= > Accuracy[%] Out[4]= > -6 Due to the manual section 3.1.6: When you do calculations with arbitrary-precision > numbers, as > discussed in the previous section, Mathematica > always keeps track of > the precision of your results, and gives only > those digits which are > known to be correct, given the precision of your > input. When you do > calculations with machine-precision numbers, > however, Mathematica > always gives you a machine[CapitalEth]precision result, > whether or not all the > digits in the result can, in fact, be determined > to be correct on the > basis of your input. Because I started with arbitrary-precision numbers > Mathematica should display > only those digits that are correct, that is none. No, 26 digits are correct Here is the number: -0.8273960599468213681 Here is the same number computed by Mathematica with 26 correct digits: -4.78339168666055402578083604864320577443814[Times]10^32 It looks like I have been using some wrong definition of correct.:-) You just proved that Precision is useless as a measure how good your numerical result is. > (check Precision instead > of Accuracy to see > this). You appear to be showing output in InputForm. If you > use OutputForm or > StandardForm only 26 digits will be shown. 32 > Out[3]= -4.7833916866605540257808360 10 InputForm is showing more because it exposes bad > digits as well as > good ones. > To relax a bit, set a new input cell to > StandardForm and type > 77617.000000000000000000000000000000000 Convert it to InputForm. You get > > 77616.999999999999999999999999999999999999999999952771`37.9031 Convert back to StandardForm > > 77616.99999999999999999999999999999999999999999976637`37.9031 Again to InputForm > > 77616.99999999999999999999999999999999999999999963735`37.9031 Back to StandardForm > > 77616.99999999999999999999999999999999999999999951376`37.9031 See what you can get if you have enough patience > or a small program. PK Agreed, it's not very pretty. I am uncertain as to > whether this > indicates a bug in StandardForm or elsewhere in the > underlying numerics > code, and will defer to our numerics experts on that > issue. My guess is > it is a bug if only because it violates the spirit > of IEEE arithmetic > wherein floats that have integer values should be > representable as such > (or something to that effect). I will point out, > however, that the two > numbers in question are equal to the specified > precision. Also it > appears to be improved in our development kernel. > Daniel Lichtblau > Wolfram Research __________________________________________________ Do you Yahoo!? http://faith.yahoo.com ==== > You are entitled to your opinion. For my applications > this behavior IS useless. > I agree that Mathematica is probably useless for you. This is however not because it is useless or useless for your application, but because to use its full power you have to study it, understand it, and in particular, for numerical work, understand the model of arithmetic it uses. Lie with mathematics they are really no shortcuts that will lead you to its full power. In addition, since it is a computer program, it has certain conventions, which may not be the same as the conventions of other programs (they all have conventions) but which you have to accept to be able to use it. Now, once you have done that, you may still not like the way Mathematica does things and there are genuine experts in numerics who indeed do not like and are quite vocal about it. But they never say it is useless, because by saying that you are either displaying your own ignorance or engaging in stupid and pointless abuse. On a more serious level, there seem to be two basic approaches to numeric computation relevant to this discussion. It seems to me (though I am no expert in this sort of thing) that there are three types of situations that one may encounter. Firstly, there is the vast majority of rather simple computations for which built-in machine floating point arithmetic , which carries no guarantee of precision at all is meant for. It clearly must be sufficient for the majority of purposes, since most general purpose and even technical software available uses not other method. The reason of course is that it is by far the fastest way to do such computations (as well as being sufficient for most situations). The second type of situation is when you actually know the precision of your input and would like the program to give you some idea about the precision of the output you might expect. This is the most likely situation in empirical science and is exactly what SetPrecision is meant for. Most reasonable people would agree that Mathematica works well in this situation. There is finally one more situation, to which the only reasonable criticism that I have read in this thread appears to be directed at. That is the situation when you actually know your input exactly, but working with exact numbers is far too slow. So what you have to do is to replace your exact numbers with inexact ones padded with 0's. In Mathematica you have to take a guess at how much padding you will need, than use SetPrecision to pad the numbers, and then check the Precision of your answer. It may turn out that you did not get as much precision as you needed, in which case you have to use more zeros. Or it may be that you used more than enough, which mans that your computation could have been done faster. I learned from Leszek Sczaniecki that there is an approach due to Oliver Aberth which lets you only specify the desired precision of your answer and the program itself will choose the correct padding for your input. It woudl ertainly be nice to have this ability, but I honestly think that it would be only of marginal advantage over making your own guess. It seems to me that the checking that the Aberth mthod must require will be time consuming and wiht a bit of practice one can probably get better results as far as speed is concerned using the Mathematica approach. But this is just pure speculation and certainly it woudl be nice if such a possibility existed. Andrzej Kozlowski Toyama International University JAPAN http://sigma.tuins.ac.jp/~andrzej/ ==== >Yes, there seems to be a lot of people who have a visceral hatred for >Microsoft and Windows. They are even willing to shed blood to avoid >Windows. But why? Windows works and you don't have to become a systems >programmer. > >Furthermore, I think that Steven Wolfram uses some version of Windows. >So guess which system Mathematica will be best tuned up for? If it is true Wolfram uses Mathematica on a Windows based machine my experience is it doesn't translate to Mathematica running better on Windows. I use Mathematica on WindowsNT and on a Mac (currently Mac OS X). I have found Mathematica to be more stable on a Mac than on Windows. On more than one occassion I've seen errors I made Mathematica code to crash the entire machine under WindowsNT. I've never had this happen running things on a Mac. ==== >Yes, there seems to be a lot of people who have a visceral hatred for >>Microsoft and Windows. They are even willing to shed blood to avoid >>Windows. But why? Windows works and you don't have to become a systems >>programmer. >> >>Furthermore, I think that Steven Wolfram uses some version of Windows. >>So guess which system Mathematica will be best tuned up for? If it is true Wolfram uses Mathematica on a Windows based machine my > experience is it doesn't translate to Mathematica running better on > Windows. I use Mathematica on WindowsNT and on a Mac (currently Mac OS X). > I have found Mathematica to be more stable on a Mac than on Windows. On > more than one occassion I've seen errors I made Mathematica code to crash > the entire machine under WindowsNT. I've never had this happen running > things on a Mac. I believe you hit the nail on the head. I suspect Dr. Wolfram is running on Mac. I have the feeling WRI is a clandestine Mac holdout. STH . ==== I ran Turbo XML http://www.tibco.com/solutions/products/extensibility/turbo_xml.jsp on ToFileName[{$TopDirectory, SystemFiles, IncludeFiles, XML}, notebookml1.dtd ] It gave me an error saying: There is more than one attribute named class. My guess is this was the intention: hattons@ljosalfr:~/.Mathematica/SystemFiles/IncludeFiles/XML/NotebookML1> diff /opt/Wolfram/Mathematica/4.2/SystemFiles/IncludeFiles/XML/NotebookML1/notebo okml.dtd /home/hattons/.Mathematica/SystemFiles/IncludeFiles/XML/NotebookML1/notebook ml.dtd 91c91 Comments? STH . ==== I've posted Mathematica notebooks and packages illustrating most of the neural networks discussed in James Anderson's book An Introduction to Neural Networks to the Brainstage Research web site (www.brainstage.com). Have fun! Don Donald Doherty, Ph.D. Brainstage Research, Inc. donald.doherty@brainstage.com ==== can I have mathematica solver things like a(over)b, (in how many ways can you pick b items from a items)? I have mathematica 4. Stefan ==== Stefan You certainly can - try Binomial[a,b]. Mark Westwood > can I have mathematica solver things like a(over)b, (in how many ways > can you pick b items from a items)? I have mathematica 4. Stefan ==== There are two basic ways, the second of which has two forms. The basic ways are: 1. Using the built in function ContourPlot, e.g.: ContourPlot[x^3y + y^3 - 9, {x, -9, 9}, {y, -27, 27}, Contours -> {0}, ContourShading -> False, Axes -> True, Frame -> False, PlotPoints -> 50, AxesOrigin -> {0, 0}] alternatively you can use a Standard package: <{0,0}] or ImplicitPlot[x^3y+y^3-==9,{x,-9,9},{y,-27,27},AxesOrigin->{0,0}] The difference between these two is that the first one gives you a smoother picture but requires the equation to be solvable (by Mathematica) for y. The second will give a picture very similar to that produced by the first method. Andrzej Kozlowski Yokohama, Japan http://www.mimuw.edu.pl/~akoz/ http://platon.c.u-tokyo.ac.jp/andrzej/ > How can I plot functions like: > > (x-2)^2 + 2(y-3)^2 = 6 > > and > > x^3y + y^3 = 9 > > using Mathematica? > > > > > > ==== Since you already know the answer, the simplest way is: In[51]:= Factor[x^4 + x^3 + x^2 + x + 1, Extension -> {GoldenRatio}] Out[51]= (-(-1 - x + GoldenRatio*x - x^2))*(1 + GoldenRatio*x + x^2) Andrzej Kozlowski Yokohama, Japan http://www.mimuw.edu.pl/~akoz/ http://platon.c.u-tokyo.ac.jp/andrzej/ > Greetings MathGroup, > > My name is Steve Earth, and I am a new subscriber to this list and > also a > new user of Mathematica; so please forgive this rather simple > question... > > I would like to enter the quartic x^4 + x^3 + x^2 + x + 1 into > Mathematica > and have it be able to tell me that it factors into > > (x^2 + GoldenRatio x + 1) ( x^2 - 1/GoldenRatio x + 1) > > What instructions do I need to execute to achieve this output? > > -Steve Earth > Harker School > http://www.harker.org/ > > > > Greetings MathGroup, > > My name is Steve Earth, and I am a new subscriber to this list and > also a > new user of Mathematica; so please forgive this rather simple > question... > > I would like to enter the quartic x^4 + x^3 + x^2 + x + 1 into > Mathematica > and have it be able to tell me that it factors into > > (x^2 + GoldenRatio x + 1) ( x^2 - 1/GoldenRatio x + 1) > > What instructions do I need to execute to achieve this output? > > -Steve Earth > Harker School > http://www.harker.org/ > > > Andrzej Kozlowski Yokohama, Japan http://www.mimuw.edu.pl/~akoz/ http://platon.c.u-tokyo.ac.jp/andrzej/ ==== I confess I am not 100% sure what you mean. Would you like to do this in steps, like you would do it by hand? In[1]:= a + b/c + d*(Sqrt[e]/c) == 0; In[2]:= Thread[(#1 - a - b/c & )[%], Equal] Out[2]= (d*Sqrt[e])/c == -a - b/c In[3]:= Thread[(#1*c & )[%], Equal] Out[3]= d*Sqrt[e] == (-a - b/c)*c In[4]:= Thread[(#1^2 & )[%], Equal] Out[4]= d^2*e == (-a - b/c)^2*c^2 In[5]:= Simplify[%] Out[5]= d^2*e == (b + a*c)^2 Of course you can combine all the steps into a single function, but I think it will be fairly complicated. My own favourite way to do this sort of thing is: In[1]:= Simplify[d^2*e == (d^2*e /. AlgebraicRules[ a + b/c + d*(Sqrt[e]/c) == 0, e])] Out[1]= d^2*e == (b + a*c)^2 However, AlgebraicRules has not been documented since version 4. It should be possible to do this using PolynomialReduce but it seems to require the sort of skill only Daniel Lichtblau possesses;) Andrzej Kozlowski Toyama International University JAPAN http://sigma.tuins.ac.jp/~andrzej/ > I'm a newbie and, of course, the first thing I want to do is apparently > one of the most complicated... > > I have an expression that looks like this: > > A + B/C + D*Sqrt[E]/C = 0 > > A,B,C,D, & E are all polynomials in x > I want it to look like this > > (D^2)*E = (A*C + B)^2 > > At that point, I'll have polynomials in x on both sides. Finally, I > want the equation to be written out with terms grouped by powers of x, > but I think I can do that part :) > > I'll be very grateful to anyone who can give me some pointers. Or, at > least point me to some tutorial in the Mathematica documentation. I've > been looking over the documentation and I found Appendix A.5 in The > Mathematica Book, but that doesn't help me. I _need_ some examples. I > did find a couple of well-written posts in this newsgroup, but not > quite > close enough to what I want. > > > Troy. > > =-=-=-=-=-=-=-=-=-= > > FYI, here's the expression I'm working with. > > > denom = Sqrt[(B^2 - r^2)^2 + 4*(r^2)*(b^2)] > cnu = (2*b^2 - B^2 + r^2)/denom > snu = -2*b*Sqrt[B^2 - b^2]/denom > sif = 2*r*b/denom > cif = (r^2 - B^2)/denom > > pdr = -Cos[ds]*Sin[q]*(snu*cif + > cnu*sif) - Sin[ds]*(cnu*cif - snu*sif) > > 0 == -(B^2 - b^2)*V^2/(r^2) + (((B*V)^2)/( > r^2) - 2*w*b*V*Cos[q]*Cos[ds] + (w* > r)^2 - (w*r*pdr)^2)*(Cos[qr])^2 > > Although I said it's a polynomial in x, it's really a polynomial in b > that I'm after. > > > > ==== <{-1,5}] ImplicitPlot[x^3y + y^3 == 9, {x, -10, 10}] Meilleures salutations Florian Jaccard -----Message d'origine----- Envoy.8e : dim., 6. octobre 2002 11:34 Ë : mathgroup@smc.vnet.net Objet : Plotting ellipses and other functions How can I plot functions like: (x-2)^2 + 2(y-3)^2 = 6 and x^3y + y^3 = 9 using Mathematica? ==== >How can I plot functions like: > >(x-2)^2 + 2(y-3)^2 = 6 > >and > >x^3y + y^3 = 9 > >using Mathematica? > Needs[Graphics`ImplicitPlot`]; ImplicitPlot[(x - 2)^2 + 2(y - 3)^2 == 6, {x, -1, 5}, {y, 1, 5}]; ImplicitPlot[x^3 y + y^3 == 9, {x, -6, 6}, {y, -6, 6}]; Bob Hanlon ==== >Could somebody please inform me how to Round numbers to a >certain Accuracy using Mathematica 4.2. This is not as easy as it >sounds. >Every function that I have read Rounds the Display, and not the actual >number. myRound[x_, n_] := Round[10^n*x]/10.^n; Table[myRound[Random[], 3], {10}] {0.044, 0.019, 0.738, 0.298, 0.917, 0.171, 0.021, 0.314, 0.658, 0.153} Bob Hanlon Reply-To: tgarza01@prodigy.net.mx ==== You might use ImplicitPlot: In[1]:= << Graphics`ImplicitPlot` In[2]:= eqn1 = (x - 2)^2 + 2*(y - 3)^2 == 6; In[3]:= ImplicitPlot[eqn1, {x, -2, 6}]; In[4]:= eqn2 = x^3*y + y^3 == 9; In[5]:= ImplicitPlot[eqn2, {x, -8, 8}]; Tomas Garza Mexico City Original Message: ----------------- ==== You may use Binomial. It could also be useful to look at the AddOn package DiscreteMath`Combinatorica`, where you will find a wealth of interesting things related to that. In[1]:= Binomial[6,2] Out[1]= 15 Tomas Garza Mexico City Original Message: ----------------- ==== Troy, True, interactive manipulation can be difficult. However, here is one way to do what you want. We have to do the same thing to both sides of the equation. (# - D*Sqrt[K]/C)&/@(A+B/C+D*Sqrt[K]/C[Equal]0 A + B/C == -((D*Sqrt[K])/C) Together/@% (B + A*C)/C == -((D*Sqrt[K])/C) #C&/@% B + A*C == (-D)*Sqrt[K] #^2&/@% (B + A*C)^2 == D^2*K NOTES. Here is how #C&/@ (lhs ==rhs) works: #C&/@ (lhs ==rhs) --> #C&[lhs]==#C&[rhs] --> lhs C == rhs C --> ... f/@( expr) is special for for Map[f, expr] expr& is special for Function[expr] Please look up Map and Function in the Help Browser. -- Allan --------------------- Allan Hayes Mathematica Training and Consulting Leicester UK www.haystack.demon.co.uk hay@haystack.demon.co.uk Voice: +44 (0)116 271 4198 > I'm a newbie and, of course, the first thing I want to do is apparently > one of the most complicated... > > I have an expression that looks like this: > > A + B/C + D*Sqrt[E]/C = 0 > > A,B,C,D, & E are all polynomials in x > I want it to look like this > > (D^2)*E = (A*C + B)^2 > > At that point, I'll have polynomials in x on both sides. Finally, I > want the equation to be written out with terms grouped by powers of x, > but I think I can do that part :) > > I'll be very grateful to anyone who can give me some pointers. Or, at > least point me to some tutorial in the Mathematica documentation. I've > been looking over the documentation and I found Appendix A.5 in The > Mathematica Book, but that doesn't help me. I _need_ some examples. I > did find a couple of well-written posts in this newsgroup, but not quite > close enough to what I want. > > > Troy. > > =-=-=-=-=-=-=-=-=-= > > FYI, here's the expression I'm working with. > > > denom = Sqrt[(B^2 - r^2)^2 + 4*(r^2)*(b^2)] > cnu = (2*b^2 - B^2 + r^2)/denom > snu = -2*b*Sqrt[B^2 - b^2]/denom > sif = 2*r*b/denom > cif = (r^2 - B^2)/denom > > pdr = -Cos[ds]*Sin[q]*(snu*cif + > cnu*sif) - Sin[ds]*(cnu*cif - snu*sif) > > 0 == -(B^2 - b^2)*V^2/(r^2) + (((B*V)^2)/( > r^2) - 2*w*b*V*Cos[q]*Cos[ds] + (w* > r)^2 - (w*r*pdr)^2)*(Cos[qr])^2 > > Although I said it's a polynomial in x, it's really a polynomial in b > that I'm after. > ==== > Troy, > True, interactive manipulation can be difficult. > However, here is one way to do what you want. > We have to do the same thing to both sides of the equation. (# - D*Sqrt[K]/C)&/@(A+B/C+D*Sqrt[K]/C[Equal]0 A + B/C == -((D*Sqrt[K])/C) I think I have to apologize for the lack of clarity in my original post. I had tried to word it carefully, but I deceived myself. I should have said: I have an expression that can be put into this form: A + B/C + D*Sqrt[K]/C = 0 A,B,C,D, & K are all polynomials in x I need to get it into that form and, in the end, I want it to look like this (D^2)*K = (A*C + B)^2 I think I gave the impression that I have polynomials A,B,C,D, & K at my fingertips. I don't. The expression I have is given at the end of this message. I'm still trying to digest the respones I've garned so far. In the meantime, I decided to post this clarification. > I'm a newbie and, of course, the first thing I want to do is >> apparently one of the most complicated... >> >> I have an expression that looks like this: >> >> A + B/C + D*Sqrt[K]/C = 0 >> >> A,B,C,D, & K are all polynomials in x >> I want it to look like this >> >> (D^2)*K = (A*C + B)^2 >> >> At that point, I'll have polynomials in x on both sides. Finally, I >> want the equation to be written out with terms grouped by powers of x, >> but I think I can do that part :) >> >> I'll be very grateful to anyone who can give me some pointers. Or, >> at least point me to some tutorial in the Mathematica documentation. >> I've been looking over the documentation and I found Appendix A.5 in >> The Mathematica Book, but that doesn't help me. I _need_ some >> examples. I did find a couple of well-written posts in this newsgroup, >> but not quite close enough to what I want. >> >> >> Troy. >> >> =-=-=-=-=-=-=-=-=-= >> >> FYI, here's the expression I'm working with. >> >> >> denom = Sqrt[(B^2 - r^2)^2 + 4*(r^2)*(b^2)] >> cnu = (2*b^2 - B^2 + r^2)/denom >> snu = -2*b*Sqrt[B^2 - b^2]/denom >> sif = 2*r*b/denom >> cif = (r^2 - B^2)/denom >> >> pdr = -Cos[ds]*Sin[q]*(snu*cif + >> cnu*sif) - Sin[ds]*(cnu*cif - snu*sif) >> >> 0 == -(B^2 - b^2)*V^2/(r^2) + (((B*V)^2)/( >> r^2) - 2*w*b*V*Cos[q]*Cos[ds] + (w* >> r)^2 - (w*r*pdr)^2)*(Cos[qr])^2 >> >> Although I said it's a polynomial in x, it's really a polynomial in >> b that I'm after. >> ==== Only on second reading I noticed the part about a,b,c,d being polynomials in x. Both methods will still work if first perform the same operation as below and finally use the replacement rule %/.{a->p[x],b->q[x],c->r[x],d->u[x],e->v[x]}, where p[x] etc are the given polynomials. Of course collecting of terms can be done with Collect[%,x]. Andrzej Kozlowski Yokohama, Japan http://www.mimuw.edu.pl/~akoz/ http://platon.c.u-tokyo.ac.jp/andrzej/ > I confess I am not 100% sure what you mean. Would you like to do this > in steps, like you would do it by hand? > > In[1]:= > a + b/c + d*(Sqrt[e]/c) == 0; > > In[2]:= > Thread[(#1 - a - b/c & )[%], Equal] > > Out[2]= > (d*Sqrt[e])/c == -a - b/c > > In[3]:= > Thread[(#1*c & )[%], Equal] > > Out[3]= > d*Sqrt[e] == (-a - b/c)*c > > In[4]:= > Thread[(#1^2 & )[%], Equal] > > Out[4]= > d^2*e == (-a - b/c)^2*c^2 > > In[5]:= > Simplify[%] > > Out[5]= > d^2*e == (b + a*c)^2 > > Of course you can combine all the steps into a single function, but I > think it will be fairly complicated. > > My own favourite way to do this sort of thing is: > > In[1]:= > Simplify[d^2*e == (d^2*e /. AlgebraicRules[ > a + b/c + d*(Sqrt[e]/c) == 0, e])] > > Out[1]= > d^2*e == (b + a*c)^2 > > However, AlgebraicRules has not been documented since version 4. It > should be possible to do this using PolynomialReduce but it seems to > require the sort of skill only Daniel Lichtblau possesses;) > > Andrzej Kozlowski > Toyama International University > JAPAN > http://sigma.tuins.ac.jp/~andrzej/ > > > > > >> I'm a newbie and, of course, the first thing I want to do is >> apparently >> one of the most complicated... >> >> I have an expression that looks like this: >> >> A + B/C + D*Sqrt[E]/C = 0 >> >> A,B,C,D, & E are all polynomials in x >> I want it to look like this >> >> (D^2)*E = (A*C + B)^2 >> >> At that point, I'll have polynomials in x on both sides. Finally, I >> want the equation to be written out with terms grouped by powers of x, >> but I think I can do that part :) >> >> I'll be very grateful to anyone who can give me some pointers. Or, at >> least point me to some tutorial in the Mathematica documentation. >> I've >> been looking over the documentation and I found Appendix A.5 in The >> Mathematica Book, but that doesn't help me. I _need_ some examples. I >> did find a couple of well-written posts in this newsgroup, but not >> quite >> close enough to what I want. >> >> >> Troy. >> >> =-=-=-=-=-=-=-=-=-= >> >> FYI, here's the expression I'm working with. >> >> >> denom = Sqrt[(B^2 - r^2)^2 + 4*(r^2)*(b^2)] >> cnu = (2*b^2 - B^2 + r^2)/denom >> snu = -2*b*Sqrt[B^2 - b^2]/denom >> sif = 2*r*b/denom >> cif = (r^2 - B^2)/denom >> >> pdr = -Cos[ds]*Sin[q]*(snu*cif + >> cnu*sif) - Sin[ds]*(cnu*cif - snu*sif) >> >> 0 == -(B^2 - b^2)*V^2/(r^2) + (((B*V)^2)/( >> r^2) - 2*w*b*V*Cos[q]*Cos[ds] + (w* >> r)^2 - (w*r*pdr)^2)*(Cos[qr])^2 >> >> Although I said it's a polynomial in x, it's really a polynomial in >> b >> that I'm after. >> >> >> >> > > ==== > Yes, there seems to be a lot of people who have a visceral hatred for > Microsoft and Windows. They are even willing to shed blood to avoid > Windows. But why? Windows works and you don't have to become a systems > programmer. Furthermore, I think that Steven Wolfram uses some version of Windows. So > guess which system Mathematica will be best tuned up for? I have no problems with Mathematica and Windows on my single computer. > There may be reasons for using a non-Microsoft operating system. But if > you are going to do it, make certain that they are good reasons. that of Windows XP by orders of magnitude. I recall when I first started been using Windows NT since October of 1992. (Yes, I know it hadn't been paper on the architecture of NT. In 1997 I was well on my way to being an MCSE. no stinkin' GUI' to 'have a look at the KDE project'. I took the latter route. The KDE has gone from a simple graphical desktop with a few more features than the CDE, (and a lot more glitches) to being the best desktop available. It's growth seems to be exponential. Windows seems, at best, to be linear. All of these are usability issues. There is another reason I don't like using Microsoft products. I've also been using Mozilla since 1995. (Yes, it has always been called Mozilla.) I was one of the original beta testers for the Netscape line of internet servers. When I saw what Netscape Communications were aiming for, Windows quickly lost its luster. Netscape products were designed from the ground up with portability in mind. They were striving for uniform functionality across all platforms. I also saw what Microsoft did to undermine Netscape's R & D resources. Microsoft would condescend to having not competition in their market. Where I come from, people don't put up with that. Where do I come from? I was born in Illinois. I'm obviously not of the opinion that closed source is unacceptable. I wouldn't be using Mathematica if I were. I suspect one day Mathematica will face a real open source challeng. Her name is Charolette. She is the mother of Mozilla. That will probably be years from now. WRI need to be prepared to adjust to that eventuality when it comes. > David Park > djmp@earthlink.net > http://home.earthlink.net/~djmp/ STH . ==== I feel as if I've finally had the breakthrough in intuitively understanding the Mathematica editor, or at least the basics. I want to explain to others what they really need to know about the editor to use it for basic purposes. Part of the reason I now understand the editor better is that I've since worked with LyX http://www.devel.lyx.org and XEmacs http://www.xemacs.org, and I've also become proficient with DocBook XML. When I first started using the editor, reading the Mathematica Help didn't seem to help. It seemed to tell me a whole lot more than I needed to know, and didn't tell me what I really needed to know. Now I'm in the position of being able to use it, but not being able to explain exactly what it is I've learned. Is there any documentation directed toward the beginner, which tells him or her what to do, what to expect, what quirks to be aware of, and etc.? I'm looking for something along the lines of click here to make this happen. If you see this, it may seem weird, but that's normal. Those brackedt on the left really mean... . STH . ==== >I would like to enter the quartic x^4 + x^3 + x^2 + x + 1 into Mathematica >and have it be able to tell me that it factors into > >(x^2 + GoldenRatio x + 1) ( x^2 - 1/GoldenRatio x + 1) > >What instructions do I need to execute to achieve this output? > soln = Factor[x^4 + x^3 + x^2 + x + 1, Extension -> GoldenRatio] // Simplify (x^2 - GoldenRatio*x + x + 1)*(x^2 + GoldenRatio*x + 1) soln = Simplify /@ (soln /. -GoldenRatio -> -1 - 1/GoldenRatio) (x^2 - x/GoldenRatio + 1)*(x^2 + GoldenRatio*x + 1) soln // FunctionExpand // FullSimplify x^4 + x^3 + x^2 + x + 1 Bob Hanlon ==== Stefan, I don't think I completely understand your question, but I think you are looking for the Binomial function in Mathematica. Binomial[4, 2] 6 David Park djmp@earthlink.net http://home.earthlink.net/~djmp/ Sender: steve@smc.vnet.net Approved: Steven M. Christensen , Moderator ==== Steve, You could use the Extension feature of Factor as documented in Help. expr = x^4 + x^3 + x^2 + x + 1 ans = Factor[expr, Extension -> {1/GoldenRatio}] (-(1/4))*(-2 - x + Sqrt[5]*x - 2*x^2)* (2 + x + Sqrt[5]*x + 2*x^2) You could also use... Factor[expr, Extension -> {Sqrt[5]}] It took me some effort to figure out how to manipulate the answer into your form. ans /. {x + Sqrt[5]*x -> (2*GoldenRatio)*x, -x + Sqrt[5]*x -> (2/GoldenRatio)*x} % /. (-4^(-1))*a_*b_ :> Simplify[-a/2]*Simplify[b/2] (-(1/4))*(-2 + (2*x)/GoldenRatio - 2*x^2)* (2 + 2*GoldenRatio*x + 2*x^2) (1 - x/GoldenRatio + x^2)*(1 + GoldenRatio*x + x^2) David Park djmp@earthlink.net http://home.earthlink.net/~djmp/ http://www.harker.org/ ==== David, To plot equations like that, simply use ImplicitPlot. Needs[Graphics`ImplicitPlot`] ImplicitPlot[(x - 2)^2 + 2(y - 3)^2 == 6, {x, -1, 5}, {y, 1, 5}]; ImplicitPlot[x^3*y + y^3 == 9, {x, -10, 10}, {y, -10, 10}]; You may have to fish a little to obtain the appropriate x and y ranges. Start by making them larger and then narrow down to the region that you want. I have put a new package at my web site for solving conic section problems in the plane. You can solve for complete information on any conic section and obtain a parametric representation for plotting it. The package also comes with complete Help documentation and examples. Using your first example (the second is not a conic). Needs[ConicSections`ConicSections`] eqn = (x - 2)^2 + 2(y - 3)^2 == 6; The routine ParseConic will take any quadratic equation and return the scale a, eccentricity e, a parametrization, and rotation matrix P, translation T and reflection matrix R that transforms the conic from standard position to its actual position. (In standard position the conic has its foci and verticies on the x-axis with the center at zero.) {{a, e}, curve[t_], {P, T, R}} = ParseConic[eqn] {{Sqrt[6], 1/Sqrt[2]}, {2 + Sqrt[6]*Cos[t], 3 + Sqrt[3]*Sin[t]}, {{{1, 0}, {0, 1}}, {2, 3}, {{1, 0}, {0, 1}}}} We could then plot the curve using ParametricPlot, which is more efficient and controllable. ParametricPlot[Evaluate[curve[t]], {t, -Pi, Pi}, AspectRatio -> Automatic, Frame -> True, Axes -> None, PlotLabel -> eqn]; Knowing a and e we can use the StandardConic routine to obtain all the information about the conic in standard position as a set of rules. standarddata = StandardConic[{a, e}] {conictype -> Ellipse, conicequation -> x^2/6 + y^2/3 == 1, coniccurve -> {Sqrt[6]*Cos[t], Sqrt[3]*Sin[t]}, coniccurvedomain -> {-Pi, Pi}, coniccenter -> {0, 0}, conicfocus -> {{Sqrt[3], 0}, {-Sqrt[3], 0}}, conicdirectrix -> {x == -2*Sqrt[3], x == 2*Sqrt[3]}, conicvertex -> {{Sqrt[6], 0}, {-Sqrt[6], 0}}} routine to obtain the same information for the conic in its actual position. standarddata // TransformEllipseRules[P, T, R] {conictype -> Ellipse, conicequation -> (1/6)*((-2 + x)^2 + 2*(-3 + y)^2) == 1, coniccurve -> {2 + Sqrt[6]*Cos[t], 3 + Sqrt[3]*Sin[t]}, coniccurvedomain -> {-Pi, Pi}, coniccenter -> {2, 3}, conicfocus -> {{2 + Sqrt[3], 3}, {2 - Sqrt[3], 3}}, conicdirectrix -> {2*Sqrt[3] + x == 2, x == 2*(1 + Sqrt[3])}, conicvertex -> {{2 + Sqrt[6], 3}, {2 - Sqrt[6], 3}}} David Park djmp@earthlink.net http://home.earthlink.net/~djmp/ Sender: steve@smc.vnet.net Approved: Steven M. Christensen , Moderator Reply-To: ==== Try this: A + B/C + D*(Sqrt[E]/C) == 0 (#1 - %[[1,{1, 2}]] & ) /@ % (C*#1 & ) /@ % (#1^2 & ) /@ % Simplify[%] Also, be aware that E is the natural logarithm base, reserved for that purpose. DrBob -----Original Message----- want the equation to be written out with terms grouped by powers of x, but I think I can do that part :) I'll be very grateful to anyone who can give me some pointers. Or, at least point me to some tutorial in the Mathematica documentation. I've been looking over the documentation and I found Appendix A.5 in The Mathematica Book, but that doesn't help me. I _need_ some examples. I did find a couple of well-written posts in this newsgroup, but not quite close enough to what I want. Troy. =-=-=-=-=-=-=-=-=-= FYI, here's the expression I'm working with. denom = Sqrt[(B^2 - r^2)^2 + 4*(r^2)*(b^2)] cnu = (2*b^2 - B^2 + r^2)/denom snu = -2*b*Sqrt[B^2 - b^2]/denom sif = 2*r*b/denom cif = (r^2 - B^2)/denom pdr = -Cos[ds]*Sin[q]*(snu*cif + cnu*sif) - Sin[ds]*(cnu*cif - snu*sif) 0 == -(B^2 - b^2)*V^2/(r^2) + (((B*V)^2)/( r^2) - 2*w*b*V*Cos[q]*Cos[ds] + (w* r)^2 - (w*r*pdr)^2)*(Cos[qr])^2 Although I said it's a polynomial in x, it's really a polynomial in b that I'm after. Reply-To: ==== It's a little ugly, but here's my solution: x^4 + x^3 + x^2 + x + 1 Simplify@Factor[%, Extension -> {GoldenRatio, 1/GoldenRatio}] Collect[%[[2]]/3, x]Collect[%[[3]]/3, x] % /. Sqrt[5] -> 2GoldenRatio - 1 % // Simplify Collect[#, x] & /@ % % /. {1 - GoldenRatio -> -1/GoldenRatio} FullSimplify@Expand@% (The last line is a check.) Bobby -----Original Message----- (x^2 + GoldenRatio x + 1) ( x^2 - 1/GoldenRatio x + 1) What instructions do I need to execute to achieve this output? -Steve Earth Harker School http://www.harker.org/ ==== >How can I plot functions like: > >(x-2)^2 + 2(y-3)^2 = 6 > >and > >x^3y + y^3 = 9 > >using Mathematica? Use ImplicitPlot in the package Graphics`ImplicitPlot` ==== David, <{0,0}] -----Original Message----- Sender: steve@smc.vnet.net Approved: Steven M. Christensen , Moderator Reply-To: jcd@q-e-d.org ==== I look for a way to right align entries output by MatrixForm. There is an ad hoc option, but I haven't figured out how it is supposed to work in Mathematica 4.0 (if it works at all). m={{MatrixForm[{0,0,0}],-123456789},{123456789,1}}; MatrixForm[m,TableAlignments->{Right,Top}] appears identical to MatrixForm[m] (with no SetOptions override) Strangely, Mathematica doesn't complain if you input ill options like: MatrixForm[m,TableAlignments->{WhateverIsIllegal,Right,Top,MakeSomeSpiralOfI t}] In Mathematica 2.2x, the behavior is rather strange: one option alone works, but when two are given (e.g. {Right,Top}), only the last seems to be acted upon. I didn't bother to try hard with the older version, so don't flame me if I'm wrong. At the other hand, when a TableAlignments->Right is given to TableForm, entries are right-aligned correctly. Can someone provide a way out (even a slow external module) or tell me how to use this option. I apologize for using a phony address in an attempt ==== I decided to export one of my notebooks to html/mathml http://baldur.globalsymmetry.com/proprietary/com/wri/notebooks/essential/ 66.92.149.152 baldur.globalsymmetry.com When I try to view it with Mozilla, I get a '?' in place of the imaginary number symbol. When I first load the page with Mozilla, I get an error telling me To properly display the MathML on this page you need to install the following fonts: CMSY 10, CMEX 10, Math1, Math2, Math4. For further infromation see: http://www.mozilla.org/projects/mathml/fonts I did what the instructions at the URL told me, and I still get the error. I looked through the fonts in the Mathematica font directory, and I found: Mathematica1Mono.9.bdf -wri-mathematica1mono-medium-r-normal--9-90-75-75-m-50-adobe-fontspecific Mathematica3.12.bdf -wri-mathematica3-medium-r-normal--12-120-75-75-p-70-adobe-fontspecific Mathematica3.24.bdf -wri-mathematica3-medium-r-normal--24-240-75-75-p-130-adobe-fontspecific Mathematica3.36.bdf -wri-mathematica3-medium-r-normal--36-360-75-75-p-210-adobe-fontspecific .... Mathematica7.12.bdf -wri-mathematica7-medium-r-normal--12-120-75-75-p-40-adobe-fontspecific Mathematica7.24.bdf .... -wri-mathematica6-medium-r-normal--12-120-75-75-p-30-adobe-fontspecific Mathematica6.24.bdf -wri-mathematica6-medium-r-normal--24-240-75-75-p-70-adobe-fontspecific Mathematica6.36.bdf -wri-mathematica6-medium-r-normal--15-150-75-75-p-50-adobe-fontspecific Mathematica5.12.bdf -wri-mathematica5-medium-r-normal--12-120-75-75-p-50-adobe-fontspecific Mathematica5.24.bdf -wri-mathematica5-medium-r-normal--24-240-75-75-p-110-adobe-fontspecific Mathematica5.36.bdf -wri-mathematica5-medium-r-normal--36-360-75-75-p-160-adobe-fontspecific Mathematica5.13.bdf -wri-mathematica5-medium-r-normal--13-130-75-75-p-50-adobe-fontspecific Mathematica4.18.bdf -wri-mathematica4-medium-r-normal--18-180-75-75-p-130-adobe-fontspecific Mathematica5.10.bdf .... I suspect this is what mozilla is looking for, but I don't know exactly how to tell it as much. Any ideas? STH . Reply-To: ==== I suppose my silliness is understandable, in light of all the confusion, both here and in the Browser on what significance arithmetic is, what bignums and bigfloats might be, etc. If many smart people are confused, there's a possibility --- just a possibility, mind you --- that it isn't entirely their fault. Yes? No? >>like 71 above, or -5 for Accuracy in the example that fooled me), but it's not a big deal For people who don't understand it as well as you, yes, it's a big deal. My purpose in all this is to understand the issue well enough to know how to proceed. I think that I do, now, but I doubt everybody on the list does. Daniel says this will all be clearer in the next release, and that's good! Bobby -----Original Message----- > > -1.180591620717411303424`71.0721*^21 > 71 > > 71.0721 digits of precision? I don't think so!! Either I am it altogether or you are just simply beating to death the point that in case of machine arithmetic (only!) Precision and Accuracy are purely formal and essentially meaningless. One can argue whether in this case there is any point of returning any value for Precision, or Accuracy (like 71 above, or -5 for Accuracy in the example that fooled me), but it's not a big deal and it most certainly does not make SetPrecision meaningless. On the contrary, SetPrecision is very useful and in fact it is SetPrecision that can tell you that the answer above is meaningless: In[8]:= f = SetAccuracy[333.75*b^6 + a^2*(11*a^2*b^2 - b^6 - 121*b^4 - 2) + 5.5*b^8 + a/(2*b), 50]; a=SetPrecision[77617.,$MachinePrecision]; b = SetPrecision[ 33096.,$MachinePrecision]; In[10]:= {f,Precision[f]} Out[10]= {1.19801754103509`0*^19, 0} I would say this is correct and show that SetPrecision is very useful indeed. It tells you (what of course you ought to already know in this case anyway) that machine precision will not give you a realiable answer in this case. If you know your numbers with a great deal of accuracy you can get an accurate answer: In[24]:= f = SetAccuracy[333.75*b^6 + a^2*(11*a^2*b^2 - b^6 - 121*b^4 - 2) + 5.5*b^8 + a/(2*b), 100]; a=SetPrecision[77617.,100]; b = SetPrecision[33096.,100]; In[26]:= {f, Precision[f]} Out[26]= {-0.82739605994682136814116509547981629199903311578438481991 781484167246798617832`61.2597, 61} Again you can be pretty sure that you got an accurate answer, provided of course your original setting of precision was valid. Honestly, to say that SetPrecision and SetAccuaracy are useless is one of the silliest thing I have read on this list in years. > Andrzej Kozlowski Yokohama, Japan http://www.mimuw.edu.pl/~akoz/ http://platon.c.u-tokyo.ac.jp/andrzej/ ==== How can I plot functions like: (x-2)^2 + 2(y-3)^2 = 6 and x^3y + y^3 = 9 using Mathematica? Reply-To: kuska@informatik.uni-leipzig.de ==== look what Graphics`ImplicitPlot` does. Jens How can I plot functions like: (x-2)^2 + 2(y-3)^2 = 6 and x^3y + y^3 = 9 using Mathematica? ==== David: For a start you can try the following: Using ContourPlot ContourPlot[(x - 2)^2 + 2(y - 3)^2 - 6, {x, -1, 5}, {y, -1, 5}, Contours -> {0}, ContourShading -> False, ContourSmoothing -> 5] Or using ImplicitPlot << Graphics`ImplicitPlot` ImplicitPlot[(x - 2)^2 + 2(y - 3)^2 == 6, {x, -1, 5}] ImplicitPlot[x^3y + y^3 == 9, {x, -4, 4}] Hans > How can I plot functions like: > > (x-2)^2 + 2(y-3)^2 = 6 > > and > > x^3y + y^3 = 9 > > using Mathematica? > > > > ==== Try reading help under the keyword ImplicitPlot | How can I plot functions like: | | (x-2)^2 + 2(y-3)^2 = 6 | | and | | x^3y + y^3 = 9 | | using Mathematica? | | | | ==== David Enquire within the Help Browser for the ImplicitPlot package and all will be revealed. For example, I snipped the following lines from the documentation: << Graphics`ImplicitPlot` ImplicitPlot[x^2 + 2 y^2 == 3, {x, -2, 2}] Hope this is of sufficient help to get you started - post again if you have any further questions. Mark Westwood How can I plot functions like: (x-2)^2 + 2(y-3)^2 = 6 and x^3y + y^3 = 9 using Mathematica? ==== you can. Use: < is there any way of linking Mathematica with Excel? Lu.92s ==== I am getting the same problem (though Mathematica 4.1). Has anyone any ideas? Aron. > >I am facing the problem in starting the link to math kernel from within >Excel. > >Specifically, when I get the message link failed to open when I click >Launch button in the 'Start Mathematica Link' dialog box. I have tried using >Multilink too so that I am able to access the kernel from Mathematica and >Excel simultaneously but I still get the same message. > >I am currently using Mathematica 4.0 and Excel 2002 (Excel XP). The programs >Mathematica and Excel otherwise appear to work fine. I am using Windows XP >home edition on Dell 8100 laptop. I have used the Mathematica Link for MS >Excel for Excel 2000 in my installation as there were no specific files for >Excel 2002. Given that I was able to successfully add the menus within excel >I think the addin should work fine but it does not? > >All help would be sincerely appreciated. > >Sincerely, > >Tahir Sheikh. ==== Download the 2002 file from this page. I use it for Excel 2002 Service Pack 2. http://support.wolfram.com/applicationpacks/excel_link/excelxp.html > > I am getting the same problem (though Mathematica 4.1). > Has anyone any ideas? > > Aron. > > > >I am facing the problem in starting the link to math kernel from within >Excel. > >Specifically, when I get the message link failed to open when I click >Launch button in the 'Start Mathematica Link' dialog box. I have tried > using >Multilink too so that I am able to access the kernel from Mathematica > and >Excel simultaneously but I still get the same message. > >I am currently using Mathematica 4.0 and Excel 2002 (Excel XP). The > programs >Mathematica and Excel otherwise appear to work fine. I am using Windows > XP >home edition on Dell 8100 laptop. I have used the Mathematica Link for > MS >Excel for Excel 2000 in my installation as there were no specific files > for >Excel 2002. Given that I was able to successfully add the menus within > excel >I think the addin should work fine but it does not? > >All help would be sincerely appreciated. > >Sincerely, > >Tahir Sheikh. > > > ==== I did not request any accuracy for f. I set the accuracy of the numerical components of the expression f. You cannot request the accuracy of the result of your computation in Mathematica, you can only set the accuracy of the input and later check what accuracy of the output results form it. In my last message on this topic I tried to explain this in the plainest and simplest way I could think of. There is nothing more left for me to say. I feel like Sisyphus but unlike him I can at least give up! Andrzej Kozlowski >> >> [...] >> >> I would say this is correct and show that SetPrecision is very useful >> indeed. It tells you (what of course you ought to already know in this >> case anyway) that machine precision will not give you a realiable >> answer in this case. If you know your numbers with a great deal of >> accuracy you can get an accurate answer: >> >> In[24]:= >> f = SetAccuracy[333.75*b^6 + a^2*(11*a^2*b^2 - b^6 - >> 121*b^4 - 2) + 5.5*b^8 + a/(2*b), 100]; >> a=SetPrecision[77617.,100]; b = SetPrecision[33096.,100]; >> >> >> In[26]:= >> {f, Precision[f]} >> >> Out[26]= >> {-0.82739605994682136814116509547981629199903311578438481991 >> 781484167246798617832`61.2597, 61} >> > > Congratulations! You just requested accuracy of 100 for f and got 61 ( > to convince yourself add Accuracy[f] to In[26]). If In[24] one > replaces SetAccuracy by SetPrecision the result is similar. > > PK > >> Again you can be pretty sure that you got an accurate answer, provided >> of course your original setting of precision was valid. >> >> Honestly, to say that SetPrecision and SetAccuaracy are useless is one >> of the silliest thing I have read on this list in years. >> >> > >> Andrzej Kozlowski >> Yokohama, Japan >> http://www.mimuw.edu.pl/~akoz/ >> http://platon.c.u-tokyo.ac.jp/andrzej/ > > > Andrzej Kozlowski Yokohama, Japan http://www.mimuw.edu.pl/~akoz/ http://platon.c.u-tokyo.ac.jp/andrzej/ Reply-To: ==== Here's a more intuitive method, perhaps: a + b/c + d*(Sqrt[e]/c) == 0 f = %[[1, 3]] %% /. f -> g First@Solve[%, g] f^2 == (g^2 /. %) However, it occurs to me you might want a more general method to collect a radical on one side and then square both sides. If so, here's a clumsy first attempt: expr = a + b/c + d*(Sqrt[e]/c) == 0 f = First@Cases[expr, Power[a_, Rational[b_, c_]], Infinity] power = First@Cases[f, Rational[b_, c_] -> Rational[c, b], Infinity] coefficient = First@Cases[expr, Times[a_, f] -> a, Infinity] Solve[expr /. coefficient f -> g, g][[1, 1]] g^2 == (g^2 /. %) % /. g -> coefficient f DrBob -----Original Message----- DrBob -----Original Message----- want the equation to be written out with terms grouped by powers of x, but I think I can do that part :) I'll be very grateful to anyone who can give me some pointers. Or, at least point me to some tutorial in the Mathematica documentation. I've been looking over the documentation and I found Appendix A.5 in The Mathematica Book, but that doesn't help me. I _need_ some examples. I did find a couple of well-written posts in this newsgroup, but not quite close enough to what I want. Troy. =-=-=-=-=-=-=-=-=-= FYI, here's the expression I'm working with. denom = Sqrt[(B^2 - r^2)^2 + 4*(r^2)*(b^2)] cnu = (2*b^2 - B^2 + r^2)/denom snu = -2*b*Sqrt[B^2 - b^2]/denom sif = 2*r*b/denom cif = (r^2 - B^2)/denom pdr = -Cos[ds]*Sin[q]*(snu*cif + cnu*sif) - Sin[ds]*(cnu*cif - snu*sif) 0 == -(B^2 - b^2)*V^2/(r^2) + (((B*V)^2)/( r^2) - 2*w*b*V*Cos[q]*Cos[ds] + (w* r)^2 - (w*r*pdr)^2)*(Cos[qr])^2 Although I said it's a polynomial in x, it's really a polynomial in b that I'm after. ==== On Friday, October 4, 2002, at 06:01 PM, DrBob > >[...] > > I would say this is correct and show that > SetPrecision is very useful > indeed. It tells you (what of course you ought > to already know in this > case anyway) that machine precision will not > give you a realiable > answer in this case. If you know your numbers > with a great deal of > accuracy you can get an accurate answer: > > In[24]:= > f = SetAccuracy[333.75*b^6 + a^2*(11*a^2*b^2 - > b^6 - > 121*b^4 - 2) + 5.5*b^8 + a/(2*b), 100]; > a=SetPrecision[77617.,100]; b = > SetPrecision[33096.,100]; > > > In[26]:= > {f, Precision[f]} > > Out[26]= > > {-0.82739605994682136814116509547981629199903311578438481991 > 781484167246798617832`61.2597, 61} > Congratulations! You just requested accuracy of > 100 for f and got 61 ( > to convince yourself add Accuracy[f] to In[26]). > If In[24] one > replaces SetAccuracy by SetPrecision the result is > similar. PK > [...] One has (initially) an accuracy of 100 for an > expression that contains > variables. In[25]:= Clear[a,b,f] In[26]:= f = SetAccuracy[333.75*b^6 + > a^2*(11*a^2*b^2 - b^6 - > 121*b^4 - 2) + 5.5*b^8 + a/(2*b), 100]; In[27]:= Accuracy[f] > Out[27]= 100. Now we assign values to some indeterminants in f. In[28]:= a = SetPrecision[77617.,100]; b = > SetPrecision[33096.,100]; In[29]:= {f, Precision[f], Accuracy[f]} > Out[29]= > {-0.8273960599468213681411650954798162919990331157843848199178148, 61.2599, 61.3422} The precision and accuracy has dropped. This is all > according to > standard numerical analysis regarding cancellation > error. You'll find it > in any textbook on the topic. > Assume that I want accuracy and precision of 100 for f. You advice me to make experiments to find out, what should be the initial precision and accuracy of a and b to reach the requested accuracy and precision for f. Notice, that you cannot just repeat I[26], we saw already what happens. I have to re-type I[24], I[25], I[26], I[27], I[28], and I[29] as many times as needed to get f with accuracy and precision 100. Dan, you simply advocate to do MANUAL WORK that should be done by machine. Let's suppose that in the above example I just want 60 digits not 61. Precisely, I want 60 digits and nothing or zeros afterwards. Let's see if I could use SetAccuracy. In[30]:= SetAccuracy[%, 60] Out[30]= -0.82739605994682136814116509547981629199903311578438481991781 In[31]:= % // FullForm Out[30]//FullForm= -0.827396059946821368141165095479816291999033115784384819917814841672467988` 59.9177 Oops, it did not work (as expected). Let's highlight with mouse the expression in Out[30] and copy to a new cell. Oops, we got -0.827396059946821368141165095479816291999033115784384819917814841672467988` 59.9177 again. Let's change Out[30] to a text cell and then copy. In[31]:= -0.82739605994682136814116509547981629199903311578438481991781 Out[31]= -0.82739605994682136814116509547981629199903311578438481991781 Success? Not so fast. In[32]:= % // FullForm Out[32]//FullForm= -0.8273960599468213681411650954798162919990331157843848199178099999999999986 35 08`59.2041 Dan, is there any simple way to get what I want? As I repeated already number of times, at this stage of the development of computer technology, software should do it for me (!). We both know that this is doable. Some of the textbooks that you just advised me to read describe it. As a developer of Mathematica, tell us why do you consider this to be a bad idea? Peter Kosta > As for what happens when you artificially raise > precision (or accuracy) > of machine numbers far beyond that guaranteed by > their internal > representation, that falls into to category of > garbage in, garbage out. > It is, howoever, valid to use SetPrecision to raise > precision in > (typically iterative) algorithms where significance > arithmetic might be > unduly pessimistic due to incorrect assumptions > about uncorollatedness > of numerical error. Examples of such usage have > appeared in this news > group. > Daniel Lichtblau > Wolfram Research __________________________________________________ Do you Yahoo!? http://faith.yahoo.com ==== > Are there any known issues with simpy treating the JLink.jar as a Java > extension as follows? > cp JLink.jar $JAVA_HOME/jre/lib/ext? > According to my understanding of the discussion in the Java Tutorial on > extensions, that should work: It does; starting with M4.2, J/Link 2.0 gets preinstalled and comes with a 1.4 murphee ==== Are there any known issues with simpy treating the JLink.jar as a Java extension as follows? cp JLink.jar $JAVA_HOME/jre/lib/ext? According to my understanding of the discussion in the Java Tutorial on extensions, that should work: http://java.sun.com/docs/books/tutorial/ext/basics/install.html Commants? STH . ==== >Are there any known issues with simpy treating the JLink.jar as a Java >extension as follows? >cp JLink.jar $JAVA_HOME/jre/lib/ext? > >According to my understanding of the discussion in the Java Tutorial on >extensions, that should work: >http://java.sun.com/docs/books/tutorial/ext/basics/install.html > >Commants? You should not do this. Code from the jre/lib/ext directory is trusted, so this poses a security risk from malicious applets. Leave JLink.jar where it lives in the JLink directory. If you want it to be available to all Java programs on your system, add its location to your CLASSPATH environment variable (this is not a security risk, as remote applets cannot load classes from CLASSPATH). Todd Gayley Wolfram Research ==== The key is in using the command Factor with the option Extension: In[1]:= Factor[x^4 + x^3 + x^2 + x + 1, Extension -> {GoldenRatio}] Out[1]= -((-1 - x + GoldenRatio*x - x^2)*(1 + GoldenRatio*x + x^2)) For manual verification you should keep in mind that: 1/GoldenRatio = GoldenRatio - 1 Germ.87n Buitrago ----- Original Message ----- > (x^2 + GoldenRatio x + 1) ( x^2 - 1/GoldenRatio x + 1) > > What instructions do I need to execute to achieve this output? > > -Steve Earth > Harker School > http://www.harker.org/ > ==== Actually, including 1/GoldenRatio in the extension leads to an unnecessarily complicated formula. In this case there is no real need to so, since by definition In[30]:= Unevaluated[1/GoldenRatio==GoldenRatio-1]//FullSimplify Out[30]= True If one really insists on having the answer in the form proposed in Steve's original posting one can simply do: (Collect[#1, x] & ) /@ Factor[x^4 + x^3 + x^2 + x + 1, Extension -> {GoldenRatio}] /. -1 + GoldenRatio -> 1/GoldenRatio (-(-1 + x/GoldenRatio - x^2))*(1 + GoldenRatio*x + x^2) > > In[]:=Factor[x^4 + x^3 + x^2 + x + 1, Extension -> {GoldenRatio, > 1/GoldenRatio}] > Out[]=-((-3 - 2*x + Sqrt[5]*x + GoldenRatio*x - 3*x^2)* > (3 + x + Sqrt[5]*x + GoldenRatio*x + 3*x^2))/9 > > Jens > >> >> Greetings MathGroup, >> >> My name is Steve Earth, and I am a new subscriber to this list and >> also a >> new user of Mathematica; so please forgive this rather simple >> question... >> >> I would like to enter the quartic x^4 + x^3 + x^2 + x + 1 into >> Mathematica >> and have it be able to tell me that it factors into >> >> (x^2 + GoldenRatio x + 1) ( x^2 - 1/GoldenRatio x + 1) >> >> What instructions do I need to execute to achieve this output? >> >> -Steve Earth >> Harker School >> http://www.harker.org/ > > > Andrzej Kozlowski Yokohama, Japan http://www.mimuw.edu.pl/~akoz/ http://platon.c.u-tokyo.ac.jp/andrzej/ ==== I want to apply a function to every k-th element of a long list and add the result to the k+1 element. [Actually k = 3 and I just want to multiply myList[[k]] by a constant (independent of k) and add the result to myList[[k+1]] for every value of k that's divisible by 3.] Is there a way to do this -- or in general to get at every k-th element of a list -- that's faster and more elegant than writing a brute force Do[] loop or using Mod[] operators, and that will take advantage of native List operators, but still not be too recondite? I've been thinking about multiplying a copy of myList by a mask list {0,0,1,0,0,1,..} to generate a masked copy and approaches like that. Better ways??? ==== Take[list, am, n, sa] gives elements m through n in steps of s. > I want to apply a function to every k-th element of a long list and > add the result to the k+1 element. [Actually k = 3 and I just want to multiply myList[[k]] by a > constant (independent of k) and add the result to myList[[k+1]] for > every value of k that's divisible by 3.] Is there a way to do this -- or in general to get at every k-th > element of a list -- that's faster and more elegant than writing a brute > force Do[] loop or using Mod[] operators, and that will take > advantage of native List operators, but still not be too recondite? I've been thinking about multiplying a copy of myList by a mask list > {0,0,1,0,0,1,..} to generate a masked copy and approaches like that. > Better ways??? ==== lst= Range[14] {1,2,3,4,5,6,7,8,9,10,11,12,13,14} A list of positions in lst ( for your purpose Range[1, Length[lst], 3] will do) ps= {3,5,10}; The following applies h to each ps element in lst and adds the result to the following element (lst[[ps+1]]=h/@lst[[ps]]+lst[[ps+1]];lst) {1,2,3,4+h[3],5,6+h[5],7,8,9,10,11+h[10],12,13,14} -- Allan --------------------- Allan Hayes Mathematica Training and Consulting Leicester UK www.haystack.demon.co.uk hay@haystack.demon.co.uk Voice: +44 (0)116 271 4198 > I want to apply a function to every k-th element of a long list and > add the result to the k+1 element. > > [Actually k = 3 and I just want to multiply myList[[k]] by a > constant (independent of k) and add the result to myList[[k+1]] for > every value of k that's divisible by 3.] > > Is there a way to do this -- or in general to get at every k-th > element of a list -- that's faster and more elegant than writing a brute > force Do[] loop or using Mod[] operators, and that will take > advantage of native List operators, but still not be too recondite? > > I've been thinking about multiplying a copy of myList by a mask list > {0,0,1,0,0,1,..} to generate a masked copy and approaches like that. > Better ways??? > ==== Consider the following approach, whish uses the MapAt command, that is Map with 'mapping-position' control. dummyFun={#,trueFun@#}& list={a,b,c,d,e,f,g,h,i,j,k,l,m,n,o,p} spec=Partition[Range[3,Length[list],3],1] MapAt[dummyFun,list,spec] %//Flatten Hope that is what you want, Borut | I want to apply a function to every k-th element of a long list and | add the result to the k+1 element. | | [Actually k = 3 and I just want to multiply myList[[k]] by a | constant (independent of k) and add the result to myList[[k+1]] for | every value of k that's divisible by 3.] | | Is there a way to do this -- or in general to get at every k-th | element of a list -- that's faster and more elegant than writing a brute | force Do[] loop or using Mod[] operators, and that will take | advantage of native List operators, but still not be too recondite? | | I've been thinking about multiplying a copy of myList by a mask list | {0,0,1,0,0,1,..} to generate a masked copy and approaches like that. | Better ways??? | Reply-To: kuska@informatik.uni-leipzig.de ==== something like: With[{k=3}, Flatten[ {#[[2]] + c*#[[1]], #[[3]]} & /@ Partition[lst, k, k, {1, 1}], 1] ] ?? Jens I want to apply a function to every k-th element of a long list and > add the result to the k+1 element. [Actually k = 3 and I just want to multiply myList[[k]] by a > constant (independent of k) and add the result to myList[[k+1]] for > every value of k that's divisible by 3.] Is there a way to do this -- or in general to get at every k-th > element of a list -- that's faster and more elegant than writing a brute > force Do[] loop or using Mod[] operators, and that will take > advantage of native List operators, but still not be too recondite? I've been thinking about multiplying a copy of myList by a mask list > {0,0,1,0,0,1,..} to generate a masked copy and approaches like that. > Better ways??? Reply-To: ==== Daniel, >>The precision/accuracy tracking mechanism will generally let you know, in some fashion, that you have no trustworthy digits. But it is up to the user to check that sort of thing. In this case Mathematica did NOT let us know, in any fashion, that we had no trustworthy digits. Precision and Accuracy outputs were completely misleading. (16 and -5 respectively.) Even Andrzej Kozlowski, who's adept in Mathematica, thought that would be meaningful, and never came up with a better way to check (other than using infinite precision for numbers that probably aren't known that exactly). Peter Kosta demonstrated that he could get a completely erroneous answer with Infinite precision. I blame the problem primarily, and I don't think there's any way to make the answer meaningful. That's not Mathematica's fault at all, and users need to be aware of that old maxim: garbage in, garbage out. comes up with a 22-digit result, it doesn't take much sophistication to realize the answer can't have 16-digit precision. Here's an even more extreme result: f = SetAccuracy[333.75*b^6 + a^2*(11*a^2*b^2 - b^6 - 121*b^4 - 2) + 5.5*b^8 + a/(2*b), 50]; a = 77617.; b = 33096.; f Precision[f] -1.180591620717411303424`71.0721*^21 71 71.0721 digits of precision? I don't think so!! We can do the following instead: x = Interval[333.75]; y = Interval[5.5]; a = Interval[77617.]; b = Interval[33096.]; x*b^6 + a^2*(11*a^2*b^2 - b^6 - 121*b^4 - 2) + y*b^8 + a/(2*b) Interval[{-4.486248158726164*^22, 4.2501298345826815*^22}] and that looks like the right answer, finally!! I like that method. However, that doesn't change the fact that Accuracy, Precision, and SetAccuracy appear to be completely useless. I haven't seen an example in which they did what anyone (but you) thought they should do. Bobby -----Original Message----- > you're not aware there's a problem, it lets you go on your merry way, > working with noise. Bobby Mathematica is not a mind reader. But the evaluation sequence, while complicated, is reasonably well documented. If you perform machine arithmetic, or for that matter significance arithmetic, and there is massive cancellation error, no use of SetAccuracy after the fact will fix it. The precision/accuracy tracking mechanism will generally let you know, in some fashion, that you have no trustworthy digits. But it is up to the user to check that sort of thing. It is not obvious to me what sort of error the software might notice to report. If you have a concise example of input, and expected output, I can look further. I've not seen anything in this thread that struck me as a failure of the software to warn the user, but maybe I missed something. Daniel ==== GentleBeings I have a straightforward implementation of successive approximations but I cannot seem to froce the code to find the correct solution when I have trig or exponentials involved. Can the assembled wisdom point to straghtforward fixes I know FindRoot works the object is to teach programming and successive approx, tho. kenf Below is the code Clear[f, g, gi, lim, r, rr, fr, gir, a, b, c, d, conv]; Plot[{x * ((x + 3)), 10*Sin[x]}, {x, 0.01, 2.4}, PlotStyle -> {{RGBColor[1, 0, 0], Thickness[ .006]}, {RGBColor[0, 0, 1], Thickness[ .006]}} ]; rr = FindRoot[x * ((x + 3)) == 10*Sin[x], {x, 2, 0.01, 2.4}]; f[a_] := a * ((a + 3)) /; a > 0; g[b_] := 10. * Sin[b] /; b > 0; gi[c_] := ArcSin[0.1*c] /; c > 0; Print[Actual root is , rr]; lim = 10; r = 2.0; conv = 10^-4; For[i = 1, i < lim, i++, { fr = f[r]; gir = gi[fr]; d = Abs[N[gir] - r]; i If[d < conv, Break[]]; r = gir; Print[The value of x = , r, found after , i, iterations,, with a tolerence , d, n] } ] Print[The value of x = , r, found after , i, iterations,, with a tolerence , d, n] Every man, woman and responsible child has an unalienable individual, civil, Constitutional and human right to obtain, own, and carry, openly or concealed, any weapon -- rifle, shotgun, handgun, machine gun, anything -- any time, any place, without asking anyone's permission. L. Neil Smith ==== > I want to apply a function to every k-th element of a long list and > add the result to the k+1 element. > > [Actually k = 3 and I just want to multiply myList[[k]] by a > constant (independent of k) and add the result to myList[[k+1]] for > every value of k that's divisible by 3.] > > Is there a way to do this -- or in general to get at every k-th > element of a list -- that's faster and more elegant than writing a > brute > force Do[] loop or using Mod[] operators, and that will take > advantage of native List operators, but still not be too recondite? > > I've been thinking about multiplying a copy of myList by a mask > list > {0,0,1,0,0,1,..} to generate a masked copy and approaches like that. > Better ways??? > > Here is a generalization of what you've asked: f[l_List,c_,k_Integer,p_Integer]:=Flatten[Block[{r=#[[k]] c},Join[Take[#,k-1],{r,#[[k+1]]+ r},Drop[#,k+1]]]&/@Partition[l,p]]/;(Mod[Length[l],p]==0&&k (Expand[#1] == Expand[#2] & )] {(1 + (1/2)*(1 - Sqrt[5])*x + x^2)* (1 + (1/2)*(1 + Sqrt[5])*x + x^2)} Andrzej Kozlowski Toyama International University JAPAN http://sigma.tuins.ac.jp/~andrzej/ > Steve > The notebook given after NOTEBOOK below contains functions for > factoring and > partial fractioning. > Here is an application to your problem: the first stage avoids our > needing > to know anything about the answer. > > fc=FactorR[x^4+x^3+x^2+x+1,x] > > (1 - (1/2)*(-1 - Sqrt[5])*x + x^2)* > (1 - (1/2)*(-1 + Sqrt[5])*x + x^2) > > Now we need to get rid of Sqrt[5] in terms of GoldenRatio. > This is rather messy: > > Simplify/@(fc/. Sqrt[5][Rule]2 GoldenRatio-1) > > (1 + x - GoldenRatio*x + x^2)*(1 + GoldenRatio*x + x^2) > > Simplify/@(%/.-GoldenRatio[Rule] 1/GoldenRatio -1) > > (1 + x/GoldenRatio + x^2)*(1 + GoldenRatio*x + x^2) > > > Another example > > PartialFractionsR[(1 + x)x/(1 - 3*x + x^2), x] > > 1 - (2*(-1 + 4*x))/((3 + Sqrt[5] - 2*x)*(-3 + 2*x)) + > (2*(-1 + 4*x))/((-3 + 2*x)*(-3 + Sqrt[5] + 2*x)) > > Simplify[%] > > (x*(1 + x))/(1 - 3*x + x^2) > > NOTEBOOK: to make a notebook from the following, copy from the next > line to > the line preceding XXX and paste into a new Mathematica notebook. > > Notebook[{ > > Cell[CellGroupData[{ > Cell[Factors and PartialFractions, Subtitle], > > Cell[Allan Hayes, 16 August 2001, Text], > > Cell[< > Here are some functions for factoring and expressing in partial > fractions over the reals and over the complex numbers. >, Text], > > Cell[BoxData[ > (Quit)], Input], > > Cell[BoxData[{ > (Off[General::spell1, General::spell]), n, > ((FactorC::usage = polynomial in x with complex coefficients, gives its factorization > over the complex numbers.n > The output may include Root objects which may be evaluated with > ToRadicals or N.>;)n), n, > ((FactorR::usage = polynomial in x with real coefficients, gives its factorization over > the reals.n > The output may include Root objects which may be evaluated with > ToRadicals or N.>;)n), n, > ((PartialFractionsC::usage = where ratl is a rational in x with complex coefficients, gives its > factorization over the complex numbers.n > The output may include Root objects which may be evaluated with > ToRadicals or N.>;)n), n, > ((PartialFractionsR::usage = where ratl is a rational in x with real coefficients, gives its > factorization over the real numbers.n > The output may include Root objects which may be evaluated with > ToRadicals or N.>;)), n, > (On[General::spell1, General::spell])}], Input, > InitializationCell->True], > > Cell[TextData[{ > FactorC[p_, x_] := , > StyleBox[(*over complex numbers*), > FontFamily->Arial, > FontWeight->Plain], > nTimes @@ Cases[Roots[p == 0, x, n Cubics -> False], u_ == > v_ -> x - v]n nFactorR[p_, x_] := , > StyleBox[(*over reals, coefficients must be real*), > FontFamily->Arial, > FontWeight->Plain], > n (Times @@ Join[Cases[#1, u_ == v_ /; Im[v] == 0 :> n x > - v], Cases[#1, u_ == v_ /; Im[v] > 0 :> n x^2 - x*2*Re[v] + > Abs[v]^2]] & )[n Roots[p == 0, x, Cubics -> False]] > }], Input, > InitializationCell->True], > > Cell[TextData[{ > PartialFractionsC[p_, x_] := , > StyleBox[(*over complex numbers*), > FontFamily->Arial, > FontWeight->Plain], > n(#+Apart[#2/FactorC[#3,x]])&@@Flatten[{PolynomialReduce[#,#2], > #2}]&[Numerator[#],Denominator[#]]&[Together[p]]n n > PartialFractionsR[p_, x_] := , > StyleBox[(*over reals, coefficients must be real*), > FontFamily->Arial, > FontWeight->Plain], > n(#+Apart[#2/FactorR[#3,x]])&@@Flatten[{PolynomialReduce[#,#2], > #2}]&[Numerator[#],Denominator[#]]&[Together[p]] > }], Input, > InitializationCell->True], > > Cell[CellGroupData[{ > > Cell[PROGRAMMING NOTES, Subsubsection], > > Cell[TextData[{ > The option , > StyleBox[Cubics->False, > FontFamily->Courier], > is used to keep the roots of cubics in , > StyleBox[Root[....], > FontFamily->Courier], > form. This is better for computation.n, > StyleBox[Re[v], > FontFamily->Courier], > and , > StyleBox[Abs[v]^2, > FontFamily->Courier], > are used rather than , > StyleBox[v+Conjugate[v] , > FontFamily->Courier], > and , > StyleBox[v*Conjugate[v], > FontFamily->Courier], > to prevent , > StyleBox[Apart, > FontFamily->Courier], > from factorising , > StyleBox[x^2 - x*2*Re[v] + Abs[v]^2], > FontFamily->Courier], > back to complex form. > }], Text] > }, Closed]], > > Cell[CellGroupData[{ > > Cell[EXAMPLES, Subsubsection], > > Cell[pol = Expand[(x - 1)*(x + 1)^2*(x^2 + x + 1)^2*(x^2 + 4)]; , > Input], > > Cell[CellGroupData[{ > > Cell[f1 = FactorC[pol, x], Input], > > Cell[BoxData[ > ((((-1) + x)) (((-2) [ImaginaryI] + > x)) ((2 [ImaginaryI] + > x)) ((1 + x))^2 (((((-1)))^(1/3) + x))^2 > (((-(((-1)))^(2/3)) + x))^2)], Output] > }, Open ]], > > Cell[CellGroupData[{ > > Cell[f2 = FactorR[pol, x], Input], > > Cell[BoxData[ > ((((-1) + x)) ((1 + x))^2 ((4 + > x^2)) ((1 + x + x^2))^2)], Output] > }, Open ]], > > Cell[CellGroupData[{ > > Cell[f3 = FactorR[x^3 + x + 1, x], Input], > > Cell[BoxData[ > (((x - Root[1 + #1 + #1^3 &, 1])) ((x^2 - > 2 x Root[(-1) + 2 #1 + 8 #1^3 &, 1] + > Root[(-1) - #1^4 + #1^6 &, 2]^2)))], Output] > }, Open ]], > > Cell[< > Root objects appear because of the option Cubics->False in Roots. > We can sometimes get radical forms, but notice the complication. >, Text], > > Cell[CellGroupData[{ > > Cell[ToRadicals[f3], Input], > > Cell[BoxData[ > (((((2/(3 (((-9) + @93)))))^(1/3) - ((1/2 > (((-9) + @93))))^(1/3)/3^(2/3) + x)) ((1/3 + > 1/3 ((29/2 - (3 @93)/2))^(1/3) + > 1/3 ((1/2 ((29 + 3 @93))))^(1/3) - > 2 ((((1/2 ((9 + @93))))^(1/3)/(2 > 3^(2/3)) - > 1/(2^(2/3) ((3 ((9 + > @93))))^(1/3)))) x + x^2)))], Output] > }, Open ]], > > Cell[Inexact forms can be found, from f3 :, Text], > > Cell[CellGroupData[{ > > Cell[N[f3], Input], > > Cell[BoxData[ > (((((0.6823278038280193`)([InvisibleSpace])) + > x)) ((((1.4655712318767682`)([InvisibleSpace])) - > 0.6823278038280193` x + x^2)))], Output] > }, Open ]], > > Cell[or directly, Text], > > Cell[CellGroupData[{ > > Cell[f3 = FactorR[x^3 + x + 1//N, x], Input], > > Cell[BoxData[ > (((((0.6823278038280193`)([InvisibleSpace])) + > x)) ((((1.4655712318767682`)([InvisibleSpace])) - > 0.6823278038280193` x + x^2)))], Output] > }, Open ]], > > Cell[Partial fractions, Text], > > Cell[CellGroupData[{ > > Cell[pf1 = PartialFractionsR[(2 + x)/pol, x], Input], > > Cell[BoxData[ > (1/(60 (((-1) + x))) - 1/(10 ((1 + x))^2) - > 39/(100 ((1 + x))) + ((-54) - 31 x)/(4225 ((4 + > x^2))) + ((-1) + 3 x)/(13 ((1 + x + x^2))^2) + (44 + > 193 x)/(507 ((1 + x + x^2))))], Output] > }, Open ]], > > Cell[CellGroupData[{ > > Cell[pf2 = PartialFractionsR[(1 + x)x/(1 - 3*x + x^2), x], > Input], > > Cell[< > 1 - (2*(-1 + 4*x))/((3 + Sqrt[5] - 2*x)*(-3 + 2*x)) + > (2*(-1 + 4*x))/((-3 + 2*x)*(-3 + Sqrt[5] + 2*x)) >, Output] > }, Open ]], > > Cell[CellGroupData[{ > > Cell[BoxData[ > (Simplify[%])], Input], > > Cell[(x*(1 + x))/(1 - 3*x + x^2), Output] > }, Open ]], > > Cell[Partial fractions will often involve Root objects , Text], > > Cell[CellGroupData[{ > > Cell[pf3 = PartialFractionsR[(1 + x)/(x^3 - x + 1), x], Input], > > Cell[BoxData[ > (((1 + > Root[1 - #1 + #1^3 &, > 1]))/((((x - > Root[1 - #1 + #1^3 &, > 1])) ((Root[1 - #1 + #1^3 &, 1]^2 - > 2 Root[1 - #1 + #1^3 &, > 1] Root[(-1) - 2 #1 + 8 #1^3 &, 1] + > Root[(-1) + #1^4 + #1^6 &, 2]^2)))) + ((x + > Root[1 - #1 + #1^3 &, 1] + > x Root[1 - #1 + #1^3 &, 1] - > 2 Root[(-1) - 2 #1 + 8 #1^3 &, 1] - > Root[(-1) + #1^4 + #1^6 &, 2]^2))/(((((-x^2) + > 2 x Root[(-1) - 2 #1 + 8 #1^3 &, 1] - > Root[(-1) + #1^4 + #1^6 &, 2]^2)) ((Root[1 - > #1 + #1^3 &, 1]^2 - > 2 Root[1 - #1 + #1^3 &, > 1] Root[(-1) - 2 #1 + 8 #1^3 &, 1] + > Root[(-1) + #1^4 + #1^6 &, 2]^2)))))], > Output] > }, Open ]], > > Cell[This can in fact be put in radical form:, Text], > > Cell[CellGroupData[{ > > Cell[ToRadicals[pf3], Input], > > Cell[BoxData[ > (((1 - ((2/(3 ((9 - @69)))))^(1/3) - ((1/2 ((9 > - @69))))^(1/3)/3^(2/3)))/(((((-(1/3)) + > 1/3 ((25/2 - (3 @69)/2))^(1/3) + > 1/3 ((1/2 ((25 + 3 @69))))^(1/3) + > (((-((2/(3 ((9 - @69)))))^(1/3)) - ((1/2 ((9 - > @69))))^(1/3)/3^(2/3)))^2 - > 2 (((-((2/(3 ((9 - @69)))))^(1/3)) - > ((1/2 ((9 - @69))))^(1/3)/3^(2/3))) ((1/24 ((864 > - 96 @69))^(1/3) + ((1/2 ((9 + @69))))^(1/3)/(2 > 3^(2/3)))))) ((((2/(3 ((9 - @69)))))^(1/3) + ((1 > /2 ((9 - @69))))^(1/3)/3^(2/3) + x)))) + ((1/3 - > 1/3 ((25/2 - (3 @69)/2))^(1/3) - ((2/(3 > ((9 - @69)))))^(1/3) - ((1/2 ((9 - @69))))^(1/3)/3 > ^(2/3) - 1/3 ((1/2 ((25 + 3 @69))))^(1/3) - > 2 ((1/24 ((864 - 96 @69))^(1/3) + ((1/2 > ((9 + @69))))^(1/3)/(2 3^(2/3)))) + > x + (((-((2/(3 ((9 - @69)))))^(1/3)) - > ((1/2 ((9 - @69))))^(1/3)/3^(2/3))) x))/(((((-(1 > /3)) + 1/3 ((25/2 - (3 @69)/2))^(1/3) + > 1/3 ((1/2 ((25 + 3 @69))))^(1/3) + > (((-((2/(3 ((9 - @69)))))^(1/3)) - ((1/2 ((9 - > @69))))^(1/3)/3^(2/3)))^2 - > 2 (((-((2/(3 ((9 - @69)))))^(1/3)) - > ((1/2 ((9 - @69))))^(1/3)/3^(2/3))) ((1/24 ((864 > - 96 @69))^(1/3) + ((1/2 ((9 + @69))))^(1/3)/(2 > 3^(2/3)))))) ((1/3 - > 1/3 ((25/2 - (3 @69)/2))^(1/3) - > 1/3 ((1/2 ((25 + 3 @69))))^(1/3) + > 2 ((1/24 ((864 - 96 @69))^(1/3) + ((1/2 > ((9 + @69))))^(1/3)/(2 3^(2/3)))) x - > x^2)))))], Output] > }, Closed]], > > Cell[We could have found the inexact form directly., Text], > > Cell[CellGroupData[{ > > Cell[BoxData[ > (PartialFractionsR[((1 + x))/((x^3 - x + 1)) // N, > x])], Input], > > Cell[BoxData[ > ((-(0.07614206365252976`/(((1.324717957244746`)( > [InvisibleSpace])) + > 1.` x))) + (((0.7982664819556426`)( > [InvisibleSpace])) + 0.07614206365252976` > x)/(((0.754877666246693`)([InvisibleSpace])) - > 1.324717957244746` x + 1.` x^2))], Output] > }, Open ]] > }, Closed]] > }, Open ]] > }, > ScreenRectangle->{{0, 1024}, {0, 709}}, > AutoGeneratedPackage->None, > WindowSize->{534, 628}, > WindowMargins->{{199, Automatic}, {0, Automatic}}, > ShowCellLabel->False, > StyleDefinitions -> Default.nb > ] > > XXX > -- > Allan > > --------------------- > Allan Hayes > Mathematica Training and Consulting > Leicester UK > www.haystack.demon.co.uk > hay@haystack.demon.co.uk > Voice: +44 (0)116 271 4198 > > >> Greetings MathGroup, >> >> My name is Steve Earth, and I am a new subscriber to this list and >> also a >> new user of Mathematica; so please forgive this rather simple >> question... >> >> I would like to enter the quartic x^4 + x^3 + x^2 + x + 1 into >> Mathematica > >> and have it be able to tell me that it factors into >> >> (x^2 + GoldenRatio x + 1) ( x^2 - 1/GoldenRatio x + 1) >> >> What instructions do I need to execute to achieve this output? >> >> -Steve Earth >> Harker School >> http://www.harker.org/ >> > > > > > > > > ==== In my earlier posting I used Union and SameTest to replace two equivalent answers (arising form the symmetry of the equation) by a single one. However, the way I did, while givign the right answer, it makes little logical sense since in general applying Expand would make any factorizations the same leaving us always with just a single one. WIthout using Union at all we get: (a + b*x + x^2)*(c + d*x + x^2) /. Select[SolveAlways[x^4 + x^3 + x^2 + x + 1 == (a + b*x + x^2)*(c + d*x + x^2), x], FreeQ[#1, I] & ] {(1 + (1/2 - Sqrt[5]/2)*x + x^2)* (1 + (1/2)*(1 + Sqrt[5])*x + x^2), (1 + (1/2)*(1 - Sqrt[5])*x + x^2)* (1 + (1/2)*(1 + Sqrt[5])*x + x^2)} Since having two identical answers differing only in the way they are written out is a bit of a nuisance, a way to get rid of one of them which does not suffer from the illogicality of my first approach is: Union[(a + b*x + x^2)*(c + d*x + x^2) /. Select[SolveAlways[x^4 + x^3 + x^2 + x + 1 == (a + b*x + x^2)*(c + d*x + x^2), x], FreeQ[#1, I] & ], SameTest -> (Simplify[First[#1] == First[#2]] & )] {(1 + (1/2)*(1 - Sqrt[5])*x + x^2)* (1 + (1/2)*(1 + Sqrt[5])*x + x^2)} > There is an equivalent approach that will give the answer without > knowing it in advance. All we need to know is that any quartic can be > factored over the reals as a product of two quadratics, so: > > > Union[(a + b*x + x^2)*(c + d*x + x^2) /. > Select[SolveAlways[x^4 + x^3 + x^2 + x + 1 == > (a + b*x + x^2)*(c + d*x + x^2), x], FreeQ[#1, I] & ], > SameTest -> (Expand[#1] == Expand[#2] & )] > > > {(1 + (1/2)*(1 - Sqrt[5])*x + x^2)* > (1 + (1/2)*(1 + Sqrt[5])*x + x^2)} > > > Andrzej Kozlowski > Toyama International University > JAPAN > http://sigma.tuins.ac.jp/~andrzej/ > > > > > >> Steve >> The notebook given after NOTEBOOK below contains functions for >> factoring and >> partial fractioning. >> Here is an application to your problem: the first stage avoids our >> needing >> to know anything about the answer. >> >> fc=FactorR[x^4+x^3+x^2+x+1,x] >> >> (1 - (1/2)*(-1 - Sqrt[5])*x + x^2)* >> (1 - (1/2)*(-1 + Sqrt[5])*x + x^2) >> >> Now we need to get rid of Sqrt[5] in terms of GoldenRatio. >> This is rather messy: >> >> Simplify/@(fc/. Sqrt[5][Rule]2 GoldenRatio-1) >> >> (1 + x - GoldenRatio*x + x^2)*(1 + GoldenRatio*x + x^2) >> >> Simplify/@(%/.-GoldenRatio[Rule] 1/GoldenRatio -1) >> >> (1 + x/GoldenRatio + x^2)*(1 + GoldenRatio*x + x^2) >> >> >> Another example >> >> PartialFractionsR[(1 + x)x/(1 - 3*x + x^2), x] >> >> 1 - (2*(-1 + 4*x))/((3 + Sqrt[5] - 2*x)*(-3 + 2*x)) + >> (2*(-1 + 4*x))/((-3 + 2*x)*(-3 + Sqrt[5] + 2*x)) >> >> Simplify[%] >> >> (x*(1 + x))/(1 - 3*x + x^2) >> >> NOTEBOOK: to make a notebook from the following, copy from the next >> line to >> the line preceding XXX and paste into a new Mathematica notebook. >> >> Notebook[{ >> >> Cell[CellGroupData[{ >> Cell[Factors and PartialFractions, Subtitle], >> >> Cell[Allan Hayes, 16 August 2001, Text], >> >> Cell[< >> Here are some functions for factoring and expressing in partial >> fractions over the reals and over the complex numbers. >>, Text], >> >> Cell[BoxData[ >> (Quit)], Input], >> >> Cell[BoxData[{ >> (Off[General::spell1, General::spell]), n, >> ((FactorC::usage = > polynomial in x with complex coefficients, gives its factorization >> over the complex numbers.n >> The output may include Root objects which may be evaluated with >> ToRadicals or N.>;)n), n, >> ((FactorR::usage = > polynomial in x with real coefficients, gives its factorization over >> the reals.n >> The output may include Root objects which may be evaluated with >> ToRadicals or N.>;)n), n, >> ((PartialFractionsC::usage = > where ratl is a rational in x with complex coefficients, gives its >> factorization over the complex numbers.n >> The output may include Root objects which may be evaluated with >> ToRadicals or N.>;)n), n, >> ((PartialFractionsR::usage = > where ratl is a rational in x with real coefficients, gives its >> factorization over the real numbers.n >> The output may include Root objects which may be evaluated with >> ToRadicals or N.>;)), n, >> (On[General::spell1, General::spell])}], Input, >> InitializationCell->True], >> >> Cell[TextData[{ >> FactorC[p_, x_] := , >> StyleBox[(*over complex numbers*), >> FontFamily->Arial, >> FontWeight->Plain], >> nTimes @@ Cases[Roots[p == 0, x, n Cubics -> False], u_ == >> v_ -> x - v]n nFactorR[p_, x_] := , >> StyleBox[(*over reals, coefficients must be real*), >> FontFamily->Arial, >> FontWeight->Plain], >> n (Times @@ Join[Cases[#1, u_ == v_ /; Im[v] == 0 :> n x >> - v], Cases[#1, u_ == v_ /; Im[v] > 0 :> n x^2 - x*2*Re[v] + >> Abs[v]^2]] & )[n Roots[p == 0, x, Cubics -> False]] >> }], Input, >> InitializationCell->True], >> >> Cell[TextData[{ >> PartialFractionsC[p_, x_] := , >> StyleBox[(*over complex numbers*), >> FontFamily->Arial, >> FontWeight->Plain], >> n(#+Apart[#2/FactorC[#3,x]])&@@Flatten[{PolynomialReduce[#,#2], >> #2}]&[Numerator[#],Denominator[#]]&[Together[p]]n n >> PartialFractionsR[p_, x_] := , >> StyleBox[(*over reals, coefficients must be real*), >> FontFamily->Arial, >> FontWeight->Plain], >> n(#+Apart[#2/FactorR[#3,x]])&@@Flatten[{PolynomialReduce[#,#2], >> #2}]&[Numerator[#],Denominator[#]]&[Together[p]] >> }], Input, >> InitializationCell->True], >> >> Cell[CellGroupData[{ >> >> Cell[PROGRAMMING NOTES, Subsubsection], >> >> Cell[TextData[{ >> The option , >> StyleBox[Cubics->False, >> FontFamily->Courier], >> is used to keep the roots of cubics in , >> StyleBox[Root[....], >> FontFamily->Courier], >> form. This is better for computation.n, >> StyleBox[Re[v], >> FontFamily->Courier], >> and , >> StyleBox[Abs[v]^2, >> FontFamily->Courier], >> are used rather than , >> StyleBox[v+Conjugate[v] , >> FontFamily->Courier], >> and , >> StyleBox[v*Conjugate[v], >> FontFamily->Courier], >> to prevent , >> StyleBox[Apart, >> FontFamily->Courier], >> from factorising , >> StyleBox[x^2 - x*2*Re[v] + Abs[v]^2], >> FontFamily->Courier], >> back to complex form. >> }], Text] >> }, Closed]], >> >> Cell[CellGroupData[{ >> >> Cell[EXAMPLES, Subsubsection], >> >> Cell[pol = Expand[(x - 1)*(x + 1)^2*(x^2 + x + 1)^2*(x^2 + 4)]; , >> Input], >> >> Cell[CellGroupData[{ >> >> Cell[f1 = FactorC[pol, x], Input], >> >> Cell[BoxData[ >> ((((-1) + x)) (((-2) [ImaginaryI] + >> x)) ((2 [ImaginaryI] + >> x)) ((1 + x))^2 (((((-1)))^(1/3) + x))^2 >> (((-(((-1)))^(2/3)) + x))^2)], Output] >> }, Open ]], >> >> Cell[CellGroupData[{ >> >> Cell[f2 = FactorR[pol, x], Input], >> >> Cell[BoxData[ >> ((((-1) + x)) ((1 + x))^2 ((4 + >> x^2)) ((1 + x + x^2))^2)], Output] >> }, Open ]], >> >> Cell[CellGroupData[{ >> >> Cell[f3 = FactorR[x^3 + x + 1, x], Input], >> >> Cell[BoxData[ >> (((x - Root[1 + #1 + #1^3 &, 1])) ((x^2 - >> 2 x Root[(-1) + 2 #1 + 8 #1^3 &, 1] + >> Root[(-1) - #1^4 + #1^6 &, 2]^2)))], Output] >> }, Open ]], >> >> Cell[< >> Root objects appear because of the option Cubics->False in Roots. >> We can sometimes get radical forms, but notice the complication. >>, Text], >> >> Cell[CellGroupData[{ >> >> Cell[ToRadicals[f3], Input], >> >> Cell[BoxData[ >> (((((2/(3 (((-9) + @93)))))^(1/3) - ((1/2 >> (((-9) + @93))))^(1/3)/3^(2/3) + x)) ((1/3 + >> 1/3 ((29/2 - (3 @93)/2))^(1/3) + >> 1/3 ((1/2 ((29 + 3 @93))))^(1/3) - >> 2 ((((1/2 ((9 + @93))))^(1/3)/(2 >> 3^(2/3)) - >> 1/(2^(2/3) ((3 ((9 + >> @93))))^(1/3)))) x + x^2)))], Output] >> }, Open ]], >> >> Cell[Inexact forms can be found, from f3 :, Text], >> >> Cell[CellGroupData[{ >> >> Cell[N[f3], Input], >> >> Cell[BoxData[ >> (((((0.6823278038280193`)([InvisibleSpace])) + >> x)) ((((1.4655712318767682`)([InvisibleSpace])) - >> 0.6823278038280193` x + x^2)))], Output] >> }, Open ]], >> >> Cell[or directly, Text], >> >> Cell[CellGroupData[{ >> >> Cell[f3 = FactorR[x^3 + x + 1//N, x], Input], >> >> Cell[BoxData[ >> (((((0.6823278038280193`)([InvisibleSpace])) + >> x)) ((((1.4655712318767682`)([InvisibleSpace])) - >> 0.6823278038280193` x + x^2)))], Output] >> }, Open ]], >> >> Cell[Partial fractions, Text], >> >> Cell[CellGroupData[{ >> >> Cell[pf1 = PartialFractionsR[(2 + x)/pol, x], Input], >> >> Cell[BoxData[ >> (1/(60 (((-1) + x))) - 1/(10 ((1 + x))^2) - >> 39/(100 ((1 + x))) + ((-54) - 31 x)/(4225 ((4 + >> x^2))) + ((-1) + 3 x)/(13 ((1 + x + x^2))^2) + (44 + >> 193 x)/(507 ((1 + x + x^2))))], Output] >> }, Open ]], >> >> Cell[CellGroupData[{ >> >> Cell[pf2 = PartialFractionsR[(1 + x)x/(1 - 3*x + x^2), x], >> Input], >> >> Cell[< >> 1 - (2*(-1 + 4*x))/((3 + Sqrt[5] - 2*x)*(-3 + 2*x)) + >> (2*(-1 + 4*x))/((-3 + 2*x)*(-3 + Sqrt[5] + 2*x)) >>, Output] >> }, Open ]], >> >> Cell[CellGroupData[{ >> >> Cell[BoxData[ >> (Simplify[%])], Input], >> >> Cell[(x*(1 + x))/(1 - 3*x + x^2), Output] >> }, Open ]], >> >> Cell[Partial fractions will often involve Root objects , Text], >> >> Cell[CellGroupData[{ >> >> Cell[pf3 = PartialFractionsR[(1 + x)/(x^3 - x + 1), x], Input], >> >> Cell[BoxData[ >> (((1 + >> Root[1 - #1 + #1^3 &, >> 1]))/((((x - >> Root[1 - #1 + #1^3 &, >> 1])) ((Root[1 - #1 + #1^3 &, 1]^2 - >> 2 Root[1 - #1 + #1^3 &, >> 1] Root[(-1) - 2 #1 + 8 #1^3 &, 1] + >> Root[(-1) + #1^4 + #1^6 &, 2]^2)))) + ((x + >> Root[1 - #1 + #1^3 &, 1] + >> x Root[1 - #1 + #1^3 &, 1] - >> 2 Root[(-1) - 2 #1 + 8 #1^3 &, 1] - >> Root[(-1) + #1^4 + #1^6 &, 2]^2))/(((((-x^2) + >> 2 x Root[(-1) - 2 #1 + 8 #1^3 &, 1] - >> Root[(-1) + #1^4 + #1^6 &, 2]^2)) ((Root[1 - >> #1 + #1^3 &, 1]^2 - >> 2 Root[1 - #1 + #1^3 &, >> 1] Root[(-1) - 2 #1 + 8 #1^3 &, 1] + >> Root[(-1) + #1^4 + #1^6 &, 2]^2)))))], >> Output] >> }, Open ]], >> >> Cell[This can in fact be put in radical form:, Text], >> >> Cell[CellGroupData[{ >> >> Cell[ToRadicals[pf3], Input], >> >> Cell[BoxData[ >> (((1 - ((2/(3 ((9 - @69)))))^(1/3) - ((1/2 ((9 >> - @69))))^(1/3)/3^(2/3)))/(((((-(1/3)) + >> 1/3 ((25/2 - (3 @69)/2))^(1/3) + >> 1/3 ((1/2 ((25 + 3 @69))))^(1/3) + >> (((-((2/(3 ((9 - @69)))))^(1/3)) - ((1/2 ((9 - >> @69))))^(1/3)/3^(2/3)))^2 - >> 2 (((-((2/(3 ((9 - @69)))))^(1/3)) - >> ((1/2 ((9 - @69))))^(1/3)/3^(2/3))) ((1/24 ((864 >> - 96 @69))^(1/3) + ((1/2 ((9 + @69))))^(1/3)/(2 >> 3^(2/3)))))) ((((2/(3 ((9 - @69)))))^(1/3) + ((1 >> /2 ((9 - @69))))^(1/3)/3^(2/3) + x)))) + ((1/3 - >> 1/3 ((25/2 - (3 @69)/2))^(1/3) - ((2/(3 >> ((9 - @69)))))^(1/3) - ((1/2 ((9 - @69))))^(1/3)/3 >> ^(2/3) - 1/3 ((1/2 ((25 + 3 @69))))^(1/3) - >> 2 ((1/24 ((864 - 96 @69))^(1/3) + ((1/2 >> ((9 + @69))))^(1/3)/(2 3^(2/3)))) + >> x + (((-((2/(3 ((9 - @69)))))^(1/3)) - >> ((1/2 ((9 - @69))))^(1/3)/3^(2/3))) x))/(((((-(1 >> /3)) + 1/3 ((25/2 - (3 @69)/2))^(1/3) + >> 1/3 ((1/2 ((25 + 3 @69))))^(1/3) + >> (((-((2/(3 ((9 - @69)))))^(1/3)) - ((1/2 ((9 - >> @69))))^(1/3)/3^(2/3)))^2 - >> 2 (((-((2/(3 ((9 - @69)))))^(1/3)) - >> ((1/2 ((9 - @69))))^(1/3)/3^(2/3))) ((1/24 ((864 >> - 96 @69))^(1/3) + ((1/2 ((9 + @69))))^(1/3)/(2 >> 3^(2/3)))))) ((1/3 - >> 1/3 ((25/2 - (3 @69)/2))^(1/3) - >> 1/3 ((1/2 ((25 + 3 @69))))^(1/3) + >> 2 ((1/24 ((864 - 96 @69))^(1/3) + ((1/2 >> ((9 + @69))))^(1/3)/(2 3^(2/3)))) x - >> x^2)))))], Output] >> }, Closed]], >> >> Cell[We could have found the inexact form directly., Text], >> >> Cell[CellGroupData[{ >> >> Cell[BoxData[ >> (PartialFractionsR[((1 + x))/((x^3 - x + 1)) // N, >> x])], Input], >> >> Cell[BoxData[ >> ((-(0.07614206365252976`/(((1.324717957244746`)( >> [InvisibleSpace])) + >> 1.` x))) + (((0.7982664819556426`)( >> [InvisibleSpace])) + 0.07614206365252976` >> x)/(((0.754877666246693`)([InvisibleSpace])) - >> 1.324717957244746` x + 1.` x^2))], Output] >> }, Open ]] >> }, Closed]] >> }, Open ]] >> }, >> ScreenRectangle->{{0, 1024}, {0, 709}}, >> AutoGeneratedPackage->None, >> WindowSize->{534, 628}, >> WindowMargins->{{199, Automatic}, {0, Automatic}}, >> ShowCellLabel->False, >> StyleDefinitions -> Default.nb >> ] >> >> XXX >> -- >> Allan >> >> --------------------- >> Allan Hayes >> Mathematica Training and Consulting >> Leicester UK >> www.haystack.demon.co.uk >> hay@haystack.demon.co.uk >> Voice: +44 (0)116 271 4198 >> >> > Greetings MathGroup, > > My name is Steve Earth, and I am a new subscriber to this list and > also a > new user of Mathematica; so please forgive this rather simple > question... > > I would like to enter the quartic x^4 + x^3 + x^2 + x + 1 into > Mathematica >> > and have it be able to tell me that it factors into > > (x^2 + GoldenRatio x + 1) ( x^2 - 1/GoldenRatio x + 1) > > What instructions do I need to execute to achieve this output? > > -Steve Earth > Harker School > http://www.harker.org/ > >> >> >> >> >> >> >> >> > > Andrzej Kozlowski Yokohama, Japan http://www.mimuw.edu.pl/~akoz/ http://platon.c.u-tokyo.ac.jp/andrzej/ Reply-To: ==== Let's be realistic. If you want 60 digits of precision, too bad! -- in the real world. There's nothing we can measure that closely. Drug concentrations in clinical trials are generally measured within 15%, for instance. Even machine precision is more than can be realistically expected in any application I can think of. Even getting a satellite to Jupiter probably involves more error in the final result than machine precision. (If not, it's because we rely on ongoing corrections and natural factors that put the satellite where it should be, such as gravity drawing it toward each rendezvous -- not on that kind of precision in propulsion or guidance.) So... unless all numerics in a problem have a theoretical origin, and could be represented in Mathematica as Infinite precision expressions... all this talk of higher-precision computation seems futile. The realistic question is this: given that I have confidence, say, in 6 digits of precision for the inputs of an expression, how many digits can I trust in the end result? Giving inputs MORE precision than they deserve isn't the best way to answer that question. Here are two methods of answering it in Mathematica. One uses bignums and the other uses Intervals. Repetitious trial and error are NOT required either way. BIGNUMS: ClearAll[a, b, f, x, y] f = x*b^6 + a^2*(11*a^2*b^2 - b^6 - 121*b^4 - 2) + y*b^8 + a/(2*b); a = SetPrecision[77617, 6]; b = SetPrecision[33096, 6]; x = SetPrecision[33375/100, 6]; y = SetPrecision[55/10, 6]; InputForm[f] -4.1396`-12.5121*^19 Several previous solutions have set the precision or accuracy of f before giving a and b (and possibly x and y) values. That results in making the exponents imprecise along with all coefficients (not just x and y), which may or may not be what we want. INTERVALS: I'll first digress to figure out what Interval is equivalent to 6-digit precision. You might not actually do this if you like the Interval method, but you have to decide SOMEHOW what Interval width to use. nums = {77617, 33096, 33375/100, 55/10}; (Interval[SetPrecision[#1, 6]] & ) /@ nums /. Interval[{a_, b_}] :> 2*((b - a)/(b + a)); InputForm[%] {3.2209438653903`0.207*^-6, 3.7768914672467`0.2761*^-6, 2.9260299625467`0.1652*^-6, 2.7743252840909`0.1421*^-6} For these numbers, # + Interval[{-1,1}]*#/630000& gets us very close, so I'll use that. The second method therefore is: g = #1 + Interval[{-1, 1}]*(#1/630000) & ; a = g[77617]; b = g[33096]; x = g[333.75]; y = g[5.5]; f = x*b^6 + a^2*(11*a^2*b^2 - b^6 - 121*b^4 - 2) + y*b^8 + a/(2*b) (Min[f] + Max[f])/2 Interval[{-2.136361928005054*^32, 2.1363651195928296*^32}] 1.5957938878075177*^26 Either method shows the answer has no trustworthy digits, but I think the second result is far easier to interpret. Here's another example, using the Sin function, whose derivative is Cos, whose magnitude is bounded by one. The precision of the Sin of a number should be GREATER than the precision of the number itself, especially when Cos is small. a = SetPrecision[Pi/2, 6]; InputForm[Sin[a]] 0.9999999999990905052982256654`11.6078 a = g[N[Pi/2]]; Sin[a] - 1 Interval[{-3.1084024243455137*^-12, 0}] I hope this was worthwhile to someone. DrBob -----Original Message----- >[...] > > I would say this is correct and show that > SetPrecision is very useful > indeed. It tells you (what of course you ought > to already know in this > case anyway) that machine precision will not > give you a realiable > answer in this case. If you know your numbers > with a great deal of > accuracy you can get an accurate answer: > > In[24]:= > f = SetAccuracy[333.75*b^6 + a^2*(11*a^2*b^2 - > b^6 - > 121*b^4 - 2) + 5.5*b^8 + a/(2*b), 100]; > a=SetPrecision[77617.,100]; b = > SetPrecision[33096.,100]; > > > In[26]:= > {f, Precision[f]} > > Out[26]= > > {-0.82739605994682136814116509547981629199903311578438481991 > 781484167246798617832`61.2597, 61} > Congratulations! You just requested accuracy of > 100 for f and got 61 ( > to convince yourself add Accuracy[f] to In[26]). > If In[24] one > replaces SetAccuracy by SetPrecision the result is > similar. PK > [...] One has (initially) an accuracy of 100 for an > expression that contains > variables. In[25]:= Clear[a,b,f] In[26]:= f = SetAccuracy[333.75*b^6 + > a^2*(11*a^2*b^2 - b^6 - > 121*b^4 - 2) + 5.5*b^8 + a/(2*b), 100]; In[27]:= Accuracy[f] > Out[27]= 100. Now we assign values to some indeterminants in f. In[28]:= a = SetPrecision[77617.,100]; b = > SetPrecision[33096.,100]; In[29]:= {f, Precision[f], Accuracy[f]} > Out[29]= > {-0.8273960599468213681411650954798162919990331157843848199178148, 61.2599, 61.3422} The precision and accuracy has dropped. This is all > according to > standard numerical analysis regarding cancellation > error. You'll find it > in any textbook on the topic. > Assume that I want accuracy and precision of 100 for f. You advice me to make experiments to find out, what should be the initial precision and accuracy of a and b to reach the requested accuracy and precision for f. Notice, that you cannot just repeat I[26], we saw already what happens. I have to re-type I[24], I[25], I[26], I[27], I[28], and I[29] as many times as needed to get f with accuracy and precision 100. Dan, you simply advocate to do MANUAL WORK that should be done by machine. Let's suppose that in the above example I just want 60 digits not 61. Precisely, I want 60 digits and nothing or zeros afterwards. Let's see if I could use SetAccuracy. In[30]:= SetAccuracy[%, 60] Out[30]= -0.82739605994682136814116509547981629199903311578438481991781 In[31]:= % // FullForm Out[30]//FullForm= -0.827396059946821368141165095479816291999033115784384819917814841672467 988` 59.9177 Oops, it did not work (as expected). Let's highlight with mouse the expression in Out[30] and copy to a new cell. Oops, we got -0.827396059946821368141165095479816291999033115784384819917814841672467 988` 59.9177 again. Let's change Out[30] to a text cell and then copy. In[31]:= -0.82739605994682136814116509547981629199903311578438481991781 Out[31]= -0.82739605994682136814116509547981629199903311578438481991781 Success? Not so fast. In[32]:= % // FullForm Out[32]//FullForm= -0.827396059946821368141165095479816291999033115784384819917809999999999 998635 08`59.2041 Dan, is there any simple way to get what I want? As I repeated already number of times, at this stage of the development of computer technology, software should do it for me (!). We both know that this is doable. Some of the textbooks that you just advised me to read describe it. As a developer of Mathematica, tell us why do you consider this to be a bad idea? Peter Kosta > As for what happens when you artificially raise > precision (or accuracy) > of machine numbers far beyond that guaranteed by > their internal > representation, that falls into to category of > garbage in, garbage out. > It is, howoever, valid to use SetPrecision to raise > precision in > (typically iterative) algorithms where significance > arithmetic might be > unduly pessimistic due to incorrect assumptions > about uncorollatedness > of numerical error. Examples of such usage have > appeared in this news > group. > Daniel Lichtblau > Wolfram Research __________________________________________________ Do you Yahoo!? http://faith.yahoo.com ==== Bobby, The example that I gave in my previous posting does not make my point, or it makes it in rather a hidden way, but it does show something interesting about computation with bigfloats. I'll explain what I mean by this and then give an example that does make the point directly. First, the previous example: Sin[#1]*10^#1*Log[1+10^(-#1)]&[15.9] -0.336629 Sin[#1]*10^#1*Log[1+10^(-#1)]&[SetPrecision[15.9,20]] Precision[%] -0.190858581374189370 17.7558 Sin[#1]*10^#1*Log[1+10^(-#1)]&[SetPrecision[15.9,7]] Precision[%] -0.19086 4.70309 (**) It looks as if the internal computations must be to a higher precision than 4 and that they start at SetPrecision[15.9,7]//FullForm 15.9000000000000003553`7 With MaxError = 10^-Accuracy[sp]//FullForm 1.590000000000001`*^-6 Roughly speaking, not more than the first seven digits are asserted to be correct. Now, the new example (taken from Stan Wagon, Programming Tips, Mathematica in Education and Research Volume 7, Number 2, 1988 p50) Clear[x] ser= Normal[Series[Cos[x],{x,0,200}]]; x= 75.0; ser//FullForm -2.7019882604300525`*^15 Probably not reliable. Set precision to 20: x= SetPrecision[75.0, 20]; (a=ser)//FullForm -16928.799183047`1.4688 MaxError= 10^-Accuracy[a] 575.263 Not good enough. Raise the precision to x= SetAccuracy[75.0,40]; (a=ser)//FullForm 1.0807905977573169155`7.7627 MaxError = 10^-Accuracy[a] !(1.3999657487996298`*^-6) That is 1.3999657487996298 10^-6 Good enough Allan --------------------- Allan Hayes Mathematica Training and Consulting Leicester UK www.haystack.demon.co.uk hay@haystack.demon.co.uk Voice: +44 (0)116 271 4198 ==== Bobby, You note some important points but there are situations in which raising accuracy is a practical necessity. 1) It may be that a calculation (perhaps describing a real world phenonenon) necessitates our raising the precision of inuts so as to work to high precision internally even when the output is not very sensitive to changes in input. Compare the following graphs pts = Table[({#1, Sin[#1]*10^#1*Log[1 + 10^(-#1)]} & )[x], {x, 15., 20., 0.1}]; ListPlot[pts, PlotJoined -> True] pts = Table[({#1, Sin[#1]*10^#1*Log[1 + 10^(-#1)]} & )[ SetAccuracy[x, 20]], {x, 15., 20., 0.1}]; ListPlot[pts, PlotJoined -> True] 2) Another situation in which high precision numbers are needed is when a computation with exact numbers would be slow or might use up all the memory available. We then need to replace the exact numbers with high precision bigfloats which causes the number of digits used internally to be restricted and the maximum error in the output to be reported. The N function will raise the precision of the exact numbers to try and reach the requested precision. There are also concerns with, for example, how the replacing bigfloats are related to the original numbers. How important this depends on the use to which the calculation is being put and how sensitive the output is to changes in inputs - in the plotting above it is not important. -- Allan --------------------- Allan Hayes Mathematica Training and Consulting Leicester UK www.haystack.demon.co.uk hay@haystack.demon.co.uk Voice: +44 (0)116 271 4198 > Let's be realistic. If you want 60 digits of precision, too bad! -- in > the real world. There's nothing we can measure that closely. Drug > concentrations in clinical trials are generally measured within 15%, for > instance. Even machine precision is more than can be realistically > expected in any application I can think of. Even getting a satellite to > Jupiter probably involves more error in the final result than machine > precision. (If not, it's because we rely on ongoing corrections and > natural factors that put the satellite where it should be, such as > gravity drawing it toward each rendezvous -- not on that kind of > precision in propulsion or guidance.) > > So... unless all numerics in a problem have a theoretical origin, and > could be represented in Mathematica as Infinite precision expressions... > all this talk of higher-precision computation seems futile. > > The realistic question is this: given that I have confidence, say, in 6 > digits of precision for the inputs of an expression, how many digits can > I trust in the end result? Giving inputs MORE precision than they > deserve isn't the best way to answer that question. Here are two > methods of answering it in Mathematica. One uses bignums and the > other uses Intervals. > > Repetitious trial and error are NOT required either way. > > BIGNUMS: > > ClearAll[a, b, f, x, y] > f = x*b^6 + a^2*(11*a^2*b^2 - b^6 - 121*b^4 - 2) + y*b^8 + a/(2*b); > a = SetPrecision[77617, 6]; > b = SetPrecision[33096, 6]; > x = SetPrecision[33375/100, 6]; > y = SetPrecision[55/10, 6]; > InputForm[f] > > -4.1396`-12.5121*^19 > > Several previous solutions have set the precision or accuracy of f > before giving a and b (and possibly x and y) values. That results in > making the exponents imprecise along with all coefficients (not just x > and y), which may or may not be what we want. > > INTERVALS: > > I'll first digress to figure out what Interval is equivalent to 6-digit > precision. You might not actually do this if you like the Interval > method, but you have to decide SOMEHOW what Interval width to use. > > nums = {77617, 33096, 33375/100, 55/10}; > (Interval[SetPrecision[#1, 6]] & ) /@ nums /. Interval[{a_, b_}] :> > 2*((b - a)/(b + a)); > InputForm[%] > > {3.2209438653903`0.207*^-6, > 3.7768914672467`0.2761*^-6, > 2.9260299625467`0.1652*^-6, > 2.7743252840909`0.1421*^-6} > > For these numbers, # + Interval[{-1,1}]*#/630000& gets us very close, so > I'll use that. The second method therefore is: > > g = #1 + Interval[{-1, 1}]*(#1/630000) & ; > a = g[77617]; > b = g[33096]; > x = g[333.75]; > y = g[5.5]; > f = x*b^6 + a^2*(11*a^2*b^2 - b^6 - 121*b^4 - 2) + y*b^8 + a/(2*b) > (Min[f] + Max[f])/2 > > Interval[{-2.136361928005054*^32, 2.1363651195928296*^32}] > > 1.5957938878075177*^26 > > Either method shows the answer has no trustworthy digits, but I think > the second result is far easier to interpret. > > Here's another example, using the Sin function, whose derivative is Cos, > whose magnitude is bounded by one. The precision of the Sin of a number > should be GREATER than the precision of the number itself, especially > when Cos is small. > > a = SetPrecision[Pi/2, 6]; > InputForm[Sin[a]] > > 0.9999999999990905052982256654`11.6078 > > a = g[N[Pi/2]]; > Sin[a] - 1 > > Interval[{-3.1084024243455137*^-12, 0}] > > I hope this was worthwhile to someone. > > DrBob > > -----Original Message----- > > > > > On Friday, October 4, 2002, at 06:01 PM, DrBob > >[...] > > I would say this is correct and show that > SetPrecision is very useful > indeed. It tells you (what of course you ought > to already know in this > case anyway) that machine precision will not > give you a realiable > answer in this case. If you know your numbers > with a great deal of > accuracy you can get an accurate answer: > > In[24]:= > f = SetAccuracy[333.75*b^6 + a^2*(11*a^2*b^2 - > b^6 - > 121*b^4 - 2) + 5.5*b^8 + a/(2*b), 100]; > a=SetPrecision[77617.,100]; b = > SetPrecision[33096.,100]; > > > In[26]:= > {f, Precision[f]} > > Out[26]= > > > {-0.82739605994682136814116509547981629199903311578438481991 > 781484167246798617832`61.2597, 61} > > > Congratulations! You just requested accuracy of > 100 for f and got 61 ( > to convince yourself add Accuracy[f] to In[26]). > If In[24] one > replaces SetAccuracy by SetPrecision the result is > similar. > > PK > [...] > > One has (initially) an accuracy of 100 for an > expression that contains > variables. > > In[25]:= Clear[a,b,f] > > In[26]:= f = SetAccuracy[333.75*b^6 + > a^2*(11*a^2*b^2 - b^6 - > 121*b^4 - 2) + 5.5*b^8 + a/(2*b), 100]; > > In[27]:= Accuracy[f] > Out[27]= 100. > > Now we assign values to some indeterminants in f. > > In[28]:= a = SetPrecision[77617.,100]; b = > SetPrecision[33096.,100]; > > In[29]:= {f, Precision[f], Accuracy[f]} > Out[29]= > > {-0.8273960599468213681411650954798162919990331157843848199178148, > > 61.2599, 61.3422} > > The precision and accuracy has dropped. This is all > according to > standard numerical analysis regarding cancellation > error. You'll find it > in any textbook on the topic. > > > Assume that I want accuracy and precision of 100 for > f. You advice me to make experiments to find out, what > should be the initial precision and accuracy of a and > b to reach the requested accuracy and precision for f. > Notice, that you cannot just repeat I[26], we saw > already what happens. I have to re-type I[24], I[25], > I[26], I[27], I[28], and I[29] as many times as needed > to get f with accuracy and precision 100. > > Dan, you simply advocate to do MANUAL WORK that should > be done by machine. > > Let's suppose that in the above example I just want 60 > digits not 61. Precisely, I want 60 digits and nothing > or zeros afterwards. Let's see if I could use > SetAccuracy. > > In[30]:= > SetAccuracy[%, 60] > > Out[30]= > -0.82739605994682136814116509547981629199903311578438481991781 > > In[31]:= > % // FullForm > > Out[30]//FullForm= > -0.827396059946821368141165095479816291999033115784384819917814841672467 > 988` > 59.9177 > > Oops, it did not work (as expected). Let's highlight > with mouse the expression in Out[30] and copy to a new > cell. Oops, we got > -0.827396059946821368141165095479816291999033115784384819917814841672467 > 988` > 59.9177 > again. Let's change Out[30] to a text cell and then > copy. > > In[31]:= > -0.82739605994682136814116509547981629199903311578438481991781 > > Out[31]= > -0.82739605994682136814116509547981629199903311578438481991781 > > Success? Not so fast. > > In[32]:= > % // FullForm > > Out[32]//FullForm= > -0.827396059946821368141165095479816291999033115784384819917809999999999 > 998635 > 08`59.2041 > > Dan, is there any simple way to get what I want? > > As I repeated already number of times, at this stage > of the development of computer technology, software > should do it for me (!). We both know that this is > doable. Some of the textbooks that you just advised me > to read describe it. As a developer of Mathematica, > tell us why do you consider this to be a bad idea? > > Peter Kosta > > As for what happens when you artificially raise > precision (or accuracy) > of machine numbers far beyond that guaranteed by > their internal > representation, that falls into to category of > garbage in, garbage out. > It is, howoever, valid to use SetPrecision to raise > precision in > (typically iterative) algorithms where significance > arithmetic might be > unduly pessimistic due to incorrect assumptions > about uncorollatedness > of numerical error. Examples of such usage have > appeared in this news > group. > > > Daniel Lichtblau > Wolfram Research > > > __________________________________________________ > Do you Yahoo!? > http://faith.yahoo.com > > > > Reply-To: ==== Here's a start, with a more complicated function: lst = Range[50]; f = #1^2 & ; k = 7; r = Range[k + 1, Length[lst], k]; lst[[r]] = lst[[r]] + f /@ lst[[r - 1]]; lst {1, 2, 3, 4, 5, 6, 7, 57, 9, 10, 11, 12, 13, 14, 211, 16, 17, 18, 19, 20, 21, 463, 23, 24, 25, 26, 27, 28, 813, 30, 31, 32, 33, 34, 35, 1261, 37, 38, 39, 40, 41, 42, 1807, 44, 45, 46, 47, 48, 49, 2451} or: g[lst_List, k_Integer?Positive, f_Function] := Module[ {result = lst, r = Range[k + 1, Length[lst], k]}, result[[r]] = result[[r]] + f /@ lst[[r - 1]]; result ] lst = Range[50]; lst = g[lst, 7, #1^2 & ] or: lst = Range[50]; lst + Drop[Prepend[MapIndexed[If[Mod[First@#2, k] == 0, f@#1, 0] &, lst], 0], -1] DrBob -----Original Message----- advantage of native List operators, but still not be too recondite? I've been thinking about multiplying a copy of myList by a mask list {0,0,1,0,0,1,..} to generate a masked copy and approaches like that. Better ways??? ==== I'm playing around with Mathematica just trying to see what happens if... One thing I came up with was lst1={a,b,c}; lst2={d,e,f}; Map[lst1,lst2} which resulted in the following rather unusual looking expression(?): {{a, b, c}[d], {a, b, c}[e], {a, b, c}[f]} I'm wondering if such a List represents something 'meaningful'. Any opinions? STH . ==== > I'm playing around with Mathematica just trying to see what happens if... One thing > I came up with was lst1={a,b,c}; lst2={d,e,f}; Map[lst1,lst2} which > resulted in the following rather unusual looking expression(?): > > {{a, b, c}[d], {a, b, c}[e], {a, b, c}[f]} > > I'm wondering if such a List represents something 'meaningful'. Any > opinions? Through /@ {{a, b, c}[d], {a, b, c}[e], {a, b, c}[f]} {{a[d], b[d], c[d]}, {a[e], b[e], c[e]}, {a[f], b[f], c[f]}} which might be useful. -- Allan --------------------- Allan Hayes Mathematica Training and Consulting Leicester UK www.haystack.demon.co.uk hay@haystack.demon.co.uk Voice: +44 (0)116 271 4198 ==== >I want to apply a function to every k-th element of a long list and >add the result to the k+1 element. > >[Actually k = 3 and I just want to multiply myList[[k]] by a >constant (independent of k) and add the result to myList[[k+1]] for >every value of k that's divisible by 3.] lst = Table[Random[], {20}] fact = 2 Here's a matrix method. f represents an indexed element of a matrix. f[x_,x_] := If[Mod[x,3]==0, fact, 1] f[x_, y_] := If[Mod[x,3]==0&&y[Equal]x+1, 2, 0] Create a matrix from the elements. arr = Array[f,{Length[lst], Length[lst]}]; Then matrix multiply. newlst1 = lst.arr Here's another way with the highly underrated MapIndexed. Create pairs of the nth and n-1th values. pairs = Transpose[{lst, Prepend[Drop[lst,-1],0]}] Create a function that takes the nth value (val), the n-1th value (prevval), and the index (num). multlst[{val_, prevval_}, {num_}] := Switch[Mod[num, 3], 0, fact*val, 1,val+ fact*prevval, 2, val ] Then, MapIndexed across the pairs. newlst2=MapIndexed[multlst, pairs] -------------------------------------------------------------- Omega Consulting The final answer to your Mathematica needs Spend less time searching and more time finding. http://www.wz.com/internet/Mathematica.html ==== I have two equations that I have solved for: x[n_] := 2331 + 8 n y[n_] := -3108 - 11n I want to include only solutions which are non-negative, that is x >= 0 and y >= 0. In this example we can do 2331 + 8n > = 0 and solve for n, n >= -291.375 and -3108 - 11 n >= 0 and solve for n, n <= -282.545 So we have -291.375 <= n <= -282.545. The integer solution set here is for n = {-290, -289, -288, -287, -286, -285, -284, -283}. So in this case we have 8 non-negative solutions. Given that I can supply x[n] and y[n], how do I go about finding the set n? ==== Here is a way of doing what you want (output cells are indented): <=0,y[n]>=0},n] -(2331/8) <= n <= -(3108/11) FullForm[soln] Inequality[Rational[-2331,8],LessEqual,n,LessEqual,Rational[-3108,11]] Range[Apply[Sequence,{Ceiling[soln[[1]]],Floor[soln[[-1]]]}]] {-291, -290, -289, -288, -287, -286, -285, -284, -283} -- Steve Luttrell West Malvern, UK > > I have two equations that I have solved for: > > x[n_] := 2331 + 8 n > y[n_] := -3108 - 11n > > I want to include only solutions which are non-negative, that is x >= 0 and > y >= 0. > > In this example we can do 2331 + 8n > = 0 and solve for n, n >= -291.375 > and -3108 - 11 n >= 0 and solve for n, n <= -282.545 > > So we have -291.375 <= n <= -282.545. > > The integer solution set here is for n = > {-290, -289, -288, -287, -286, -285, -284, -283}. > > So in this case we have 8 non-negative solutions. > > Given that I can supply x[n] and y[n], how do I go about finding the set n? > > > > ==== Let's first consider your original problem and take a small list as an example: mylist = {a, b, c, d, e, f, g}; As you pointed out there is a rather obvious and natural way to do it using the Do loop. In[5]:= Do[mylist[[3*i + 1]] = k*mylist[[3*i]] + mylist[[3*i + 1]], {i, 1, Length[mylist]/3}]; mylist Out[5]= {a, b, c, d + c*k, e, f, g + f*k} One can also do it using something like your second approach. Notice the need to re-set mylist which got changed by the Do loop: In[6]:= mylist={a,b,c,d,e,f,g}; In[7]:= RotateRight[Table[If[Mod[i,3]==0,k, 0],{i,1,Length[mylist]}]*mylist]+mylist Out[7]= {a,b,c,d+c k,e,f,g+f k} There is an ambiguity that appears if the length of the list is exactly divisible by 3. In that case k times the last element should be added to the next element. In this case the first approach (using Do) will produce an error message while the second will interpret the next element to mean the first element of the list. One can fix that but I shan't bother to do so and assume that the length of the list is not divisible by 3. Now let's take a large list and compare the performance of the two approaches: In[8]:= k = 2; In[9]:= mylist = Table[Random[Integer, {1, 9}], {10000}]; In[11]:= mylist1 = mylist; In[12]:= Timing[list1 = (Do[mylist[[3*i + 1]] = k*mylist[[3*i]] + mylist[[3*i + 1]], {i, 1, Length[mylist]/3}]; mylist); ] Out[12]= {0.21999999999999997*Second, Null} In[13]:= mylist = mylist1; In[14]:= Timing[list2 = RotateRight[Table[If[Mod[i, 3] == 0, k, 0], {i, 1, Length[mylist]}]*mylist] + mylist; ] Out[14]= {0.020000000000000018*Second, Null} In[15]:= list1 == list2 Out[15]= True So you were right, at least in this implementation the second approach turns out to be much faster. Note also however, that when you change the problem to the one that you originally stated the two approaches will no longer give the same answer. The reason is that your statement is ambiguous: apply a function to every k-th element of a long list and add the result to the k+1 element can either mean that you want to take the original k-th element of the original list, multiply by a constant and add to the next one, or do the same with the already altered k-th element (by the previous step of the procedure). The Do loop approach will do the latter and the other the former. Andrzej Kozlowski Yokohama, Japan http://www.mimuw.edu.pl/~akoz/ http://platon.c.u-tokyo.ac.jp/andrzej/ > I want to apply a function to every k-th element of a long list and > add the result to the k+1 element. > > [Actually k = 3 and I just want to multiply myList[[k]] by a > constant (independent of k) and add the result to myList[[k+1]] for > every value of k that's divisible by 3.] > > Is there a way to do this -- or in general to get at every k-th > element of a list -- that's faster and more elegant than writing a > brute > force Do[] loop or using Mod[] operators, and that will take > advantage of native List operators, but still not be too recondite? > > I've been thinking about multiplying a copy of myList by a mask > list > {0,0,1,0,0,1,..} to generate a masked copy and approaches like that. > Better ways??? > > > ==== Can Mathematica factor the polynomial p1=x^6+9/14*x^5+9/28*x^4+3/35*x^3+9/700*x^2+9/8750*x+3/87500; without a priori knowledge of the Extension field? ==== Carlos, Futher to my previous posting (which gave the code for the function FactorR used below), here is a complete factorisation by radicals. I also test that the product of the factors gives the original polynomial. We want to factor the polynomial p1 = x^6 + (9/14)*x^5 + (9/28)*x^4 + (3/35)*x^3 + (9/700)*x^2 + (9/8750)*x + 3/87500; in radicals. We can't expect this to be easy or even possible in terms of radicals (the general quintic is not solvable interms of radicals). But, using the function FactorR given in my posting, Re:factoring quartic over radicals, sent a few days ago (08/012/02) , we get p2 = FactorR[p1, x] (x^2 - 2*x*Root[3 + 45*#1 + 225*#1^2 + 700*#1^3 & , 1] + Root[-3 + 225*#1^2 - 5625*#1^4 + 87500*#1^6 & , 2]^2)* (x^2 - 2*x*Root[1827 + 65340*#1 + 974700*#1^2 + 7824000*#1^3 + 36360000*#1^4 + 100800000*#1^5 + 156800000*#1^6 & , 2] + Root[9 - 1350*#1^2 + 84375*#1^4 - 3056250*#1^6 - 11250000*#1^8 - 984375000*#1^10 + 7656250000*#1^12 & , 3]^2)* (x^2 - 2*x*Root[1827 + 65340*#1 + 974700*#1^2 + 7824000*#1^3 + 36360000*#1^4 + 100800000*#1^5 + 156800000*#1^6 & , 1] + Root[9 - 1350*#1^2 + 84375*#1^4 - 3056250*#1^6 - 11250000*#1^8 - 984375000*#1^10 + 7656250000*#1^12 & , 4]^2) Try to change the root objects to radical form: p3 = p2 /. r_Root :> ToRadicals[r] (3/140 + (1/140)*(13/5)^(2/3)*3^(1/3) - (1/140)*(13/5)^(1/3)*3^(2/3) - 2*(-(3/28) - (1/28)*(13/5)^(2/3)*3^(1/3) + (1/28)*(13/5)^(1/3)*3^(2/3))*x + x^2)*(x^2 - 2*x*Root[1827 + 65340*#1 + 974700*#1^2 + 7824000*#1^3 + 36360000*#1^4 + 100800000*#1^5 + 156800000*#1^6 & , 2] + Root[9 - 1350*#1 + 84375*#1^2 - 3056250*#1^3 - 11250000*#1^4 - 984375000*#1^5 + 7656250000*#1^6 & , 1])* (x^2 - 2*x*Root[1827 + 65340*#1 + 974700*#1^2 + 7824000*#1^3 + 36360000*#1^4 + 100800000*#1^5 + 156800000*#1^6 & , 1] + Root[9 - 1350*#1 + 84375*#1^2 - 3056250*#1^3 - 11250000*#1^4 - 984375000*#1^5 + 7656250000*#1^6 & , 2]) We succeeded with the first factor: f1 = p3[[1]] 3/140 + (1/140)*(13/5)^(2/3)*3^(1/3) - (1/140)*(13/5)^(1/3)*3^(2/3) - 2*(-(3/28) - (1/28)*(13/5)^(2/3)*3^(1/3) + (1/28)*(13/5)^(1/3)*3^(2/3))*x + x^2 The product of the other two factors, in a form avoiding root objects, is easily found by division: q = PolynomialQuotient[p1, f1, x] 3/3500 + ((13/5)^(1/3)*3^(2/3))/3500 + (3/175 + (1/175)*(13/5)^(1/3)*3^(2/3))*x + (9/70 - (1/140)*(13/5)^(2/3)*3^(1/3) + (1/28)*(13/5)^(1/3)*3^(2/3))*x^2 + (3/7 - (1/14)*(13/5)^(2/3)*3^(1/3) + (1/14)*(13/5)^(1/3)*3^(2/3))*x^3 + x^4 Try FactorR on this f23 = FactorR[q, x] (x^2 + ((-(1/2))*Sqrt[117/1960 + (9/784)*(13/5)^(2/3)*3^(1/3) + (39*(13/5)^(1/3)*3^(2/3))/3920] + Im[(1/2)*Sqrt[117/1960 + (9/784)*(13/5)^(2/3)*3^(1/3) + (39*(13/5)^(1/3)*3^(2/3))/3920 + (-2250 + 25*13^(2/3)*15^(1/3) - 125*13^(1/3)*15^(2/3))/8750 + (3*(7500 - 250*13^(2/3)*15^(1/3) + 250*13^(1/3)*15^(2/3))^2)/ 1225000000 - (-((2*(300 + 20*13^(1/3)*15^(2/3)))/4375) + (1/17500)*((7500 - 250*13^(2/3)*15^(1/3) + 250*13^(1/3)*15^(2/3))* ((2250 - 25*13^(2/3)*15^(1/3) + 125*13^(1/3)*15^(2/3))/4375 - (7500 - 250*13^(2/3)*15^(1/3) + 250*13^(1/3)*15^(2/3))^2/ 306250000)))/(4*Sqrt[-(117/1960) - (9/784)*(13/5)^(2/3)*3^(1/3) - (39*(13/5)^(1/3)*3^(2/3))/3920])]])^2 - 2*x*((-7500 + 250*13^(2/3)*15^(1/3) - 250*13^(1/3)*15^(2/3))/70000 + Re[(1/2)*Sqrt[117/1960 + (9/784)*(13/5)^(2/3)*3^(1/3) + (39*(13/5)^(1/3)*3^(2/3))/3920 + (-2250 + 25*13^(2/3)*15^(1/3) - 125*13^(1/3)*15^(2/3))/8750 + (3*(7500 - 250*13^(2/3)*15^(1/3) + 250*13^(1/3)*15^(2/3))^2)/ 1225000000 - (-((2*(300 + 20*13^(1/3)*15^(2/3)))/4375) + (1/17500)*((7500 - 250*13^(2/3)*15^(1/3) + 250*13^(1/3)*15^(2/3))* ((2250 - 25*13^(2/3)*15^(1/3) + 125*13^(1/3)*15^(2/3))/4375 - (7500 - 250*13^(2/3)*15^(1/3) + 250*13^(1/3)*15^(2/3))^2/ 306250000)))/(4*Sqrt[-(117/1960) - (9/784)*(13/5)^(2/3)*3^(1/3) - (39*(13/5)^(1/3)*3^(2/3))/3920])]]) + ((-7500 + 250*13^(2/3)*15^(1/3) - 250*13^(1/3)*15^(2/3))/70000 + Re[(1/2)*Sqrt[117/1960 + (9/784)*(13/5)^(2/3)*3^(1/3) + (39*(13/5)^(1/3)*3^(2/3))/3920 + (-2250 + 25*13^(2/3)*15^(1/3) - 125*13^(1/3)*15^(2/3))/8750 + (3*(7500 - 250*13^(2/3)*15^(1/3) + 250*13^(1/3)*15^(2/3))^2)/ 1225000000 - (-((2*(300 + 20*13^(1/3)*15^(2/3)))/4375) + (1/17500)*((7500 - 250*13^(2/3)*15^(1/3) + 250*13^(1/3)*15^(2/3))* ((2250 - 25*13^(2/3)*15^(1/3) + 125*13^(1/3)*15^(2/3))/4375 - (7500 - 250*13^(2/3)*15^(1/3) + 250*13^(1/3)*15^(2/3))^2/ 306250000)))/(4*Sqrt[-(117/1960) - (9/784)*(13/5)^(2/3)*3^(1/3) - (39*(13/5)^(1/3)*3^(2/3))/3920])]])^2)* (x^2 + ((1/2)*Sqrt[117/1960 + (9/784)*(13/5)^(2/3)*3^(1/3) + (39*(13/5)^(1/3)*3^(2/3))/3920] + Im[(-(1/2))*Sqrt[117/1960 + (9/784)*(13/5)^(2/3)*3^(1/3) + (39*(13/5)^(1/3)*3^(2/3))/3920 + (-2250 + 25*13^(2/3)*15^(1/3) - 125*13^(1/3)*15^(2/3))/8750 + (3*(7500 - 250*13^(2/3)*15^(1/3) + 250*13^(1/3)*15^(2/3))^2)/ 1225000000 + (-((2*(300 + 20*13^(1/3)*15^(2/3)))/4375) + (1/17500)*((7500 - 250*13^(2/3)*15^(1/3) + 250*13^(1/3)*15^(2/3))* ((2250 - 25*13^(2/3)*15^(1/3) + 125*13^(1/3)*15^(2/3))/4375 - (7500 - 250*13^(2/3)*15^(1/3) + 250*13^(1/3)*15^(2/3))^2/ 306250000)))/(4*Sqrt[-(117/1960) - (9/784)*(13/5)^(2/3)*3^(1/3) - (39*(13/5)^(1/3)*3^(2/3))/3920])]])^2 - 2*x*((-7500 + 250*13^(2/3)*15^(1/3) - 250*13^(1/3)*15^(2/3))/70000 + Re[(-(1/2))*Sqrt[117/1960 + (9/784)*(13/5)^(2/3)*3^(1/3) + (39*(13/5)^(1/3)*3^(2/3))/3920 + (-2250 + 25*13^(2/3)*15^(1/3) - 125*13^(1/3)*15^(2/3))/8750 + (3*(7500 - 250*13^(2/3)*15^(1/3) + 250*13^(1/3)*15^(2/3))^2)/ 1225000000 + (-((2*(300 + 20*13^(1/3)*15^(2/3)))/4375) + (1/17500)*((7500 - 250*13^(2/3)*15^(1/3) + 250*13^(1/3)*15^(2/3))* ((2250 - 25*13^(2/3)*15^(1/3) + 125*13^(1/3)*15^(2/3))/4375 - (7500 - 250*13^(2/3)*15^(1/3) + 250*13^(1/3)*15^(2/3))^2/ 306250000)))/(4*Sqrt[-(117/1960) - (9/784)*(13/5)^(2/3)*3^(1/3) - (39*(13/5)^(1/3)*3^(2/3))/3920])]]) + ((-7500 + 250*13^(2/3)*15^(1/3) - 250*13^(1/3)*15^(2/3))/70000 + Re[(-(1/2))*Sqrt[117/1960 + (9/784)*(13/5)^(2/3)*3^(1/3) + (39*(13/5)^(1/3)*3^(2/3))/3920 + (-2250 + 25*13^(2/3)*15^(1/3) - 125*13^(1/3)*15^(2/3))/8750 + (3*(7500 - 250*13^(2/3)*15^(1/3) + 250*13^(1/3)*15^(2/3))^2)/ 1225000000 + (-((2*(300 + 20*13^(1/3)*15^(2/3)))/4375) + (1/17500)*((7500 - 250*13^(2/3)*15^(1/3) + 250*13^(1/3)*15^(2/3))* ((2250 - 25*13^(2/3)*15^(1/3) + 125*13^(1/3)*15^(2/3))/4375 - (7500 - 250*13^(2/3)*15^(1/3) + 250*13^(1/3)*15^(2/3))^2/ 306250000)))/(4*Sqrt[-(117/1960) - (9/784)*(13/5)^(2/3)*3^(1/3) - (39*(13/5)^(1/3)*3^(2/3))/3920])]])^2) We try to get rid of the parts Re[.] and Im[.]:, f231 = f23 /. z:(_Re | _Im) :> ToRadicals[FullSimplify[z]] (((-7500 + 250*13^(2/3)*15^(1/3) - 250*13^(1/3)*15^(2/3))/70000 - (1/280)*Sqrt[3*(-390 + 13*13^(2/3)*15^(1/3) +