0
|
19476 Views
|
16 Replies
|
10 Total Likes
View groups...
Share
GROUPS:

# Why is it OK for low-precision Solve to work like this?

Posted 11 years ago
 Consider the following.Solve[SetPrecision[1234,4]/x==1]{{x->1.000}}x=1234, not x=1.  I brought this to the attention of Mathematica Support a number of times in a number of ways (I think they're sick of hearing from me), and each time they said Mathematica was functioning fine.  Mathematica Support also gave me various workarounds, e.g., use of Rationalize for use with low-precision numbers.  I'm a Home Edition user.It seems to me that Mathematica is giving the wrong answer, and it issues no warning.I showed this to a friend, and he agreed with me: something is wrong.  Would anyone else out there like to comment?We weren't just trying to find ways to make Mathematica produce questionable results.  I certainly depend on Mathematica's numeric results, and never questioned Mathematica output until I ran across a situation where it returned a negative number for a lattice constant.  This is shown in the artificial example below.  ElementData can return low-precision (measured) numbers.  In some cases, Mathematica can compute "weird" numeric outputs.  Mathematic Support said Mathematica is working fine, and it us up to me to make sure I avoid problems like this.  (Until I spotted the negative number, I didn't know there were problems like this to avoid!  That's part of the "problem": Mathematica issues no warning, as it does under other conditions in other contexts!) (* here's how I came upon this (what I call) weirdness: computing the so-called lattice constant (a) and atomic radius (r) for an atom in an FCC crystal; in the simplified example below, the lattice constant is the distance between selected atoms in a crystal *)  (* the following artificial example is for one atom at an FCC lattice point; under limited testing, until I ran into Ni, I had no reason to doubt Mathematica's output; I did not notice Ni was "weird" because of its magnitude; I noticed it was "weird" because of its sign; so, are [Pb] and Ag wrong? how do I know? *)  elemar[elem_String]:=Module[ {gramsPerAtom,atomsPerUnitCell,volumeOfUnitCell,reqn,density}, gramsPerAtom=Quantity[ElementData[elem,"AtomicWeight"],"g/mol"]/Quantity[1,"AvogadroNumber/mol"]; atomsPerUnitCell=4; volumeOfUnitCell=a^3;reqn=4r==Sqrt[2]a;density=Quantity[ElementData[elem,"Density"],"kg/m^3"];{elem,Solve[{reqn,gramsPerAtom atomsPerUnitCell/volumeOfUnitCell==density},{a,r},Reals]}];Grid[Map[elemar,{"Pb","Ag","Ni","Au"}],Alignment->Left]Pb {{a->0.0445 m,r->0.01573683247679559 m}}Ag {{a->4.088*10^-10 m,r->1.445209045881131*10^-10 m}}Ni {{a->-7.31*10^-9 m,r->-2.582782168919389*10^-9 m}}Au {{a->-1.09*10^-8 m,r->0 m}}EDIT: "Pg" above has been edited to "Pb".  Clearly the Pb values are wrong (by orders of magnitude).  And they're positive.
16 Replies
Sort By:
Posted 11 years ago
 I'll address what I can.(1) Entering input and formatting as such does have issues. I believe they are receiving some attention but, as you might guess, people hereabouts get put on multiple projects and sometimes some have to wait for others.(2) I stand corrected, it is Edit > Copy As > Input Text. Also Plain Text might work.(3) The notebook you attached is fine insofar as I, and others, can run the calculations therein without retyping. Basically any form of input that can be run without retyping is just fine.(4) I ran into a separate issue when testing it, which I sent to the relevant person. As best i can tell modulo that issue, the actual Solve precision problem seems to be fixed.(5) I am not familiar with the correspondence that went to Support. All I can say is there was a bug showing in the original post from 6 months ago. You raise the issue of how something like that can go undetected. It's because, unlike NSolve, Solve was never really intended for approximate numerical computation. Of late, with the Quantity support now in place, Solve, Integrate, and others do more frequently get used with approximate values. So unforseen problerms do appear, and we fix them when we find them.(6) The resources that go into updates are considerable. One can certainly take the position that we should put more into that. I will point out that it would happen at the expense of regular releases (which come with both new/extended functionality, and a considerable number of fixes).(7) Solving over reals is, for some classes of problem (especially multiple equations and unknowns), much harder than solving over complexes. I don't dispute what you stated about slowness in the paragraph titled "Solve and inequalities", just pointing out that there are technical reasons for the slowness. As for the issue of "nonsensical answers", if you are referiing to what you found in the original post, then yes, it's far from what we'd want. All I can say is it was not by intent. If you had in mind something else then I would need to know what it is in order to comment.
Posted 11 years ago
 Daniel Lichtblau - props for the calm response! I got angry just reading that post
Posted 11 years ago
Posted 11 years ago
 Solve is basically for infinite precision problems.   NSolve has an option, "WorkingPrecision".Solve also has that option, but I believe it is for situations involving inequalities.  The default valueis Infinity.
Posted 11 years ago
 the original source, Solve[SetPrecision[1234,4]/x==1], is of course simple to enter or to cut/paste. I tried it prior to responding. It seems fine in the development version of Mathematica.The source in the follow-up nb is not something I would care to retype. Nor would I go to the trouble of locating OCR software (or trying to code it in Mathematica). That would be very much an assymetry in our relative efforts to look into this issue. Now that I have the notebook I'll see what I can find out.As for entering input in a way that can be cut-and-pasted, I recommend as a first attempt simply copying from nb to here. Copy should convert to Mathematica's InputForm. If that fails, use Edit > Copy As > InputForm. Entering the output that way is often fine, or it can be done as an image if you prefer.
Posted 11 years ago
 I thought the source was simple enough that you could just enter it.  OCR could've been used.  Regardless, I've attached a notebook (I believe).  I also got rid of the use of my functions, to make clear what the Mma code is.  And I used a fresh kernel.  It appears to be the same basic thing as the original post.  In this case, I would like to think that 8-digit precision is "not too low."  Clearly there's something wrong with Mma, and Wolfram's agreed (based on the original post).  What I cannot fathom is why Wolfram hasn't fixed it in a timely manner nor alerted the user community.  Imagine you bought a calculator and it didn't always work?  Remember when the Intel chip had difficulty doing floating-point operations?  Am I being unreasonable?  I realize no one has  to answer me.  But can anyone explain why I think this is important and Wolfram apparently doesn't?The image below should be in the attached as a notebook.  If not, please advise, and this time tell me how you want me to upload the source.  I do know I can copy as text, and just paste that in here.  It would look like a mess, based on my experience.  I'll do that if you want me to do that. Attachments:
Posted 11 years ago
 (1) It would be helpful if this were posted in plain cut-and-pastable text. As it stands I cannot do anything to test it.(2) I do not know of any release in the past 6 months. The problems you reported might have been addressed for the next release.
Posted 11 years ago
Posted 11 years ago
 Dear robv,my name is Peter Fleck, I'm heading the Wolfram Technology Group which among other services also provides Technical Support to our valued customers.Thank you very much for sending us your issue reports to support@wolfram.com, and particularly for boiling the issue down to a simple one-line evaluation. We respond to most support requests within 3 business days, but cannot guarantee this turnaround time in all cases.You should be receiving the reply to your most recent report, the one alerting us to the Solve[SetPrecision[1234,4]/x==1] example, very soon.Thank you,Peter
Posted 11 years ago
Posted 11 years ago
 The results for In[51] and In[60]  come from an overhaul of Solve for Mathematica 8.0. Extra, unused variables are not treated the same now. Boilerplate I had:====== The difference between 7.0 and 8.0 output is that 7.0 Solve was treating equations that involved only variables as assumptions. This functionality was not precisely defined or consistently implemented and has been removed in 8.0.  Instead there is a new option MaxExtraConditions which provides a well-defined and extended version of the functionality.  ======MaxExtraConditions is given a number at least as big as the number of unused variables.  In[2]:= Solve[{r==x,12344/x==1},x, Method->"Legacy" ]                                   Out[2]= {}  In[3]:= Solve[{r==x,12344/x==1},x, MaxExtraConditions->1]                                Solve::ratnz: Solve was unable to solve the system with inexact coefficients. The answer     was obtained by solving a corresponding exact system and numericizing the result.  Out[3]= {{x -> ConditionalExpression[1234., r == 1234.]}} p.s., None of the previous repliers in this thread are in the Technical Support Group. p.p.s. The incorrect-result problem in Solve[SetPrecision[1234,4]/x==1]is different than the extra-variable issue in Solve[{r==x,12344/x==1},x]
Posted 11 years ago
 I've tried to reply to specific messages, but the Reply button doesn't always work.  I think I am replying to the whole discussion now.  I would like to address some issues.My buddy and I are aware of workarounds.  MachinePrecision and Rationalize[xxx,0] seem to work pretty well.Regarding the Solve syntax error: if that's truly the case, Mathematica should warn us.  Consider the following two examples.In[47]:= Sin[1,2]During evaluation of In[47]:= Sin::argx: Sin called with 2 arguments; 1 argument is expected. >>Out[47]= Sin[1,2]In[48]:= SetPrecision[5]During evaluation of In[48]:= SetPrecision::argr: SetPrecision called with 1 argument; 2 arguments are expected. >>Out[48]= SetPrecision[5]The much-valued online help in Mathematica does not exactly use BNF in its syntax templates, but they get the message across.  But I would like to call your attention to a feature of Solve that is not (no longer?) documented.  Parameter 3 is a list of variables you do not want to show up in the solution.    For example (NOTE: if I convert this to Code, it automatically deletes the other Code example below, so I have not converted to Code format):In[50]:= Solve[{x+y+z==1,x+2y+3z==2,x+4y+9z==3}]Out[50]= {{x->-(1/2),y->2,z->-(1/2)}}In[51]:= Solve[{x+y+z==1,x+2y+3z==2,x+4y+9z==3},{x}]Out[51]= {}In[52]:= Solve[{x+y+z==1,x+2y+3z==2,x+4y+9z==3},{x},{y,z}]Out[52]= {{x->-(1/2)}}The phenomenon at "line 51" was maddening until I wrote a function to search for solutions I wanted.  Anyway, that was a case of a parameter that is not documented.  You can tell Solve what variables not to include in its solutions.I don't remember if we noticed that the problem "seemed to go away" when we added the ",x" to the Solve.  However, that does not solve the problem in general.  Consider the following sequence, which is closer to what was going on when I discovered the "problem." In[58]:= Solve[12344/x==1] Out[58]= {{x->1.000}} In[59]:= Solve[12344/x==1,x] During evaluation of In[59]:= Solve::ratnz: Solve was unable to solve the system with inexact coefficients. The answer was obtained by solving a corresponding exact system and numericizing the result. >> Out[59]= {{x->1234.}} In[60]:= Solve[{r==x,12344/x==1},x] Out[60]= {} In[61]:= Solve[{r==x,12344/x==1},{r,x}] During evaluation of In[61]:= Solve::ratnz: Solve was unable to solve the system with inexact coefficients. The answer was obtained by solving a corresponding exact system and numericizing the result. >>Out[61]= {{r->1234.,x->1234.}}In[62]:= Solve[{r==x,12344/x==1},{r,x},Reals]Out[62]= {{r->1.000000000000000,x->1.000000000000000}}In other words, adding "Reals" results in the "problem."  We have found all sorts of ways to get weird results, but thought we boiled it down to the leading example.Finally, we (I) submitted multiple "bug" reports on this to no avail.  Mathematica Support kept telling me it was not a problem.  If Mathematica Support is now claiming there is a problem, what did I do wrong?Sorry for the sloppy exposition, but sometimes Reply doesn't work, sometimes converting to Code deletes other lines, etc.  However, publishing and then editing allowed me to "format" the Code.  Go figure...
Posted 11 years ago
 Let's say it is implicitly documented. Besides, this exact example used to work just fine in previous versions. The appropriate developers are now aware of this issue and will be looking into it.
Posted 11 years ago
 The form Solve[eqns]is undocumented, so I would not expect it to work. However, the documented formSolve[eqns, vars]works:Solve[1234`4/x == 1., x](* {{x -> 1234.}} *)Of course, that's accompanied by the warningSolve::ratnz: Solve was unable to solve the system with inexact coefficients.   The answer was obtained by solving a corresponding exact system and numericizing the result.
Posted 11 years ago
 Also I am curious that why NDSolve seems to work ok with the low precision: The last value also shows that the working precision is just around 4.
Posted 11 years ago
 No, this is definitely not OK. Thank you for bringing this issue to our attention. For now, I would suggest using a workaround along the lines of Solve[SetPrecision[{reqn, gramsPerAtom atomsPerUnitCell/volumeOfUnitCell == density}, 20], {a, r}, Reals](* Out[22]= Pb    {{a->4.951026479623944*10^-10 m,r->1.750452198788126*10^-10 m}}Ag    {{a->4.087668466298753*10^-10 m,r->1.445209045881131*10^-10 m}}Ni    {{a->3.524028160167008*10^-10 m,r->1.245932104573222*10^-10 m}}Au    {{a->4.077382180171091*10^-10 m,r->1.441572294544084*10^-10 m}} *)