@mnhoff Roundoff error, upon floating-point evaluation following numeric substitution, can vary greatly according to the particular form of the symbolic expression. Symbolic reformulation such as is done by the simplify command does not necessarily produce a rearrangement that will be prone to less roundoff, in general.

Yes, sol3 seems to suffer more at the numeric data point you mentioned. I don't see a forcing reason why this should (or should not, I'm afraid) be true for all other data inputs.

> forget(evalf):
> evalf({eval(sol1,[x1=-17.275,x2=-180.93,y1=238.59,y2=-156.49])}):
> sort(%);                                                                
          {-143.4745136, -0.01428010101, 0.01428010101, 143.4745136}

> forget(evalf):
> evalf[100]({eval(sol1,[x1=-17.275,x2=-180.93,y1=238.59,y2=-156.49])}):
> sort(evalf(%));                                                         
          {-143.4744804, -0.01427413114, 0.01427413114, 143.4744804}

> 
> forget(evalf):                                                          
> evalf(eval(sol2,[x1=-17.275,x2=-180.93,y1=238.59,y2=-156.49])):         
> sort(indets({%},identical(M)=anything));                                
  {M = -143.4584940, M = -0.01427433450, M = 0.01427433450, M = 143.4584940}

> forget(evalf):                                                   
> evalf[100](eval(sol2,[x1=-17.275,x2=-180.93,y1=238.59,y2=-156.49])):
> sort(indets({evalf(%)},identical(M)=anything));                       
  {M = -143.4744804, M = -0.01427413114, M = 0.01427413114, M = 143.4744804}

> 
> forget(evalf):  
> evalf(eval(sol3,[x1=-17.275,x2=-180.93,y1=238.59,y2=-156.49])):
> sort(indets({%},identical(M)=anything));
  {M = -295.0142973, M = -0.01425070076, M = 0.01425070076, M = 295.0142973}

> forget(evalf):                                                   
> evalf[100](eval(sol3,[x1=-17.275,x2=-180.93,y1=238.59,y2=-156.49])):
> sort(indets({evalf(%)},identical(M)=anything));
  {M = -143.4744804, M = -0.01427413114, M = 0.01427413114, M = 143.4744804}

It seems like a reasonable supposition that, "in the main", lengthier symbolic reformulations are more at risk of such float-evaluation roundoff. But it is not true in general.

ps. I've been wanting to do a followup to an earlier bloggish post relating roundoff to session dependent "order issues" (ie. about sums rather than products). But you may find this of some interest. It won't help you decide whether your `sol1` is better for uses in exported (C, Matlab, or other double-precision computational environments).

acer

What data structure do you use to handle the data? How do you add new items? Can you post any sample code?

acer

@itsme I believe that is also a problem with the usual `Explore` command. Basically, it (must) use not only a new worksheet but also a new kernel/engine. So it does not allow (lexical or other) scoping of varables to get their values from the original worksheet.The thing to be explore has to be "fully explicit" in the sense you describe.

I have been experimenting with the version of Explore which inserts into the current worksheet. That should solve the scoping problem. (Kamel B. also asked, I think...) But it's been very busy at work and so i've had very little spare time. I'll try and remember to send you a note if i get a chance to make it workable (if crude).

@itsme I believe that is also a problem with the usual `Explore` command. Basically, it (must) use not only a new worksheet but also a new kernel/engine. So it does not allow (lexical or other) scoping of varables to get their values from the original worksheet.The thing to be explore has to be "fully explicit" in the sense you describe.

I have been experimenting with the version of Explore which inserts into the current worksheet. That should solve the scoping problem. (Kamel B. also asked, I think...) But it's been very busy at work and so i've had very little spare time. I'll try and remember to send you a note if i get a chance to make it workable (if crude).

@mnhoff For the R2008b to R2009b releases of Matlab one could configure the symbolic toolbox to use Maple instead of MuPAD as the symbolic engine. See this note. One had to have a working Maple installed. Maybe you have access to such an older Matlab? (I don't know which are the Maple versions which which it would be compatible.)

@mnhoff For the R2008b to R2009b releases of Matlab one could configure the symbolic toolbox to use Maple instead of MuPAD as the symbolic engine. See this note. One had to have a working Maple installed. Maybe you have access to such an older Matlab? (I don't know which are the Maple versions which which it would be compatible.)

It does work in Maple 13, which may well be what the poster wants. Resourceful!

This is also a good time to say thanks to Markiyan Hirnyk for his significant efforts on this forum, and congratulations as today his reputation moved into the top ten.

ps. This approach doesn't work in Maple 10, but I vaguely recall some `solve` issues restricted to just that release.

acer

It does work in Maple 13, which may well be what the poster wants. Resourceful!

This is also a good time to say thanks to Markiyan Hirnyk for his significant efforts on this forum, and congratulations as today his reputation moved into the top ten.

ps. This approach doesn't work in Maple 10, but I vaguely recall some `solve` issues restricted to just that release.

acer

@Stavros Probably true. The referenced code does not properly handle the non-full-rank situations (esp. in which leading zeros in the RREF will end up in columns > m the row-dimension).

A:=Matrix(2,4,symbol=a):
A[2,1..2]:=Vector[row](2):
A;
E:=ElemDecomp(A[1..2,[1,3]]); # otherwise falls down
Einv:=seq(LinearAlgebra:-MatrixInverse(E[i],method=pseudo),i=nops([E])..1,-1);
map(normal,`.`(Einv).A);
LinearAlgebra:-LUDecomposition(A,output=R);

Easier to write one's own simple GJ, which forms (directly, without inverting) and stores the E[i]^(-1) as it goes.

@Stavros Probably true. The referenced code does not properly handle the non-full-rank situations (esp. in which leading zeros in the RREF will end up in columns > m the row-dimension).

A:=Matrix(2,4,symbol=a):
A[2,1..2]:=Vector[row](2):
A;
E:=ElemDecomp(A[1..2,[1,3]]); # otherwise falls down
Einv:=seq(LinearAlgebra:-MatrixInverse(E[i],method=pseudo),i=nops([E])..1,-1);
map(normal,`.`(Einv).A);
LinearAlgebra:-LUDecomposition(A,output=R);

Easier to write one's own simple GJ, which forms (directly, without inverting) and stores the E[i]^(-1) as it goes.

@Preben Alsholm I wanted to ensure that evalhf was not being used, even when Digits<evalhf(Digits).

I wanted to gauge behaviour when Digits is strictly adhered to in a "software float" manner, so I added a no-op to the proc which should make it non-evalhfable.

Thank you, Alejandro.

I hadn't searched for it as I feared it might non-centralized. But if that is the formula, then that would be a useful clarifying detail for the plot,computation help page.

acer

Thank you, Alejandro.

I hadn't searched for it as I feared it might non-centralized. But if that is the formula, then that would be a useful clarifying detail for the plot,computation help page.

acer

See this old post, and its comments, for some views on the computational performance section of that ncrunch review. Some of the Maple code used in the performance comparison section of that review was not optimally written. Also, the functionality sections had some errors or omissions.

What's your thesis topic? If we knew, we might be able to make better suggestions on places to look.

acer

@Christopher2222 It should not be hard, at all. The DiscreteTransforms:-ForurierTransform command will take your data Array M from the sample (time) space to a data Array P in the frequency space. (You could also window that subtask, by taking several smaller samples and averaging, etc.) You can compute the proportion of the total power using amplitudes (using `abs`, say), and that includes computing the proportion contributed by a given subrange of the frequencies. You can compute the power spectrum. You can filter out the noise by removing from P those frequencies whose amplitude is below a threshold of your discretion, eg, that which contributes only a proportion of the whole. You can even rebuild your signal, without that noise.

If you don't know how to start at all, look at the "Startup Code' in the help worksheet ?applications,SignalGeneration where you study especially the code for the routine `discrMagnitudePlot`. Note its use in,

SigData_f := DiscreteTransforms:-FourierTransform (SigData_n);
plt := discrMagnitudePlot(SigData_f, Fs);
DocumentTools:-SetProperty ('plt_freq', 'value', plt);

(It's not a great worksheet, and lacks internal comments in the code that would explain it. But together with wikipedia, you should be able to figure it out. maybe there's a more complete solution on the Aopplication Center. Sorry, I don't have the time...)