Someone constructs an ill-conditioned situation and then uses it to prove that a certain package has limitations. To me this sounds a bit like the Maple-Mma comparisons that are occasionally pilfered on the web in order to demonstrate the superiority of one over the other.

In this case, I do not believe the example says much about the quality of DirectSearch, or the other tools in Maple. If faced with an exponent of 100 in a real situation I would not even dream of using any package in a direct straightforward application. Axel gets my vote here: first bring your problem into a form amenable to a solution. If that is not possible, you are likely up to no good in the sense that the underlying problem may just be too ill-conditioned to be amenable to a solution. And granted, I see this from a physicists perspective, meaning that in my world there would usually be some underlying physics or engineering problem. If the problem is ill-conditioned, likely any solution would not be very stable in practice.

Personally, I prefer to use Maple's build-in facilities (solve and fsolve for solution, the various fiting modules for linear and nonlinear fits). Do they have limitations: you bet they do. Can one work around those: in a surprising number of cases: yes. What do I do if Maple's modules do not succeed: I go back to my problem reviewing what I am actually trying to do here. I have dabbled with DirectSearch but not seriously used it, but in the cases where I tried it its results were not out of line with the ones I found using Maple's modules. In that sense I don't think DirectSearch is any worse; and I can believe that it will solve problems that Maple's modules have a hard time with.

But finding roots involving powers of 100 using numerical methods... well, I have a hard time with that no matter what they are called and who programmed them. And it certainly does not disqualify a software routine if it can't do that without a lot of hand-holding.

Just my $0.02

M.D.

@Carl Love I tried your logarithmic plotting application. Getting it to load (in Maple 15 in this case) seemed easy enough; but I run into trouble when running logplots.mws: The first example produces no error but also no output (the linear plot works but Logplot( ) produces nothing); the 2nd example produces an error:

Error:Error creating plot:nullPlot Internal Error:unrecognised option Dag in Plot.

That unrecognized option Dag shows up repeatedly for the other examples.

If you can check this and a fix is possible with reasonable effort I would appreciate your looking into it. At the same time I understand that fixing this may require too much effort...

Many thanks,

M.D.

@tomleslie Right. Doing that... and running into the old bug in plots:-display that affects Maple 15 and Maple 17, the two I am using most of the time: Once you use plots:-setoptions() (which I have in my .mapleinit to set things like gridlines and axis label fonts), axis[n]=[mode=log] no longer works right. So I cannot get the plots to be log-log for the amplitude and log-lin for the phase.

I tried to list plots:-display with interface(verboseproc=2); but I only get this:

and I cannot list plots/display. Listing "plots" in this way just gives me the list of exports but not the actual members of "plots". My idea was to see if by any chance I could make a custom version of plots:-display that works itself around the issue with log scales. The bug was confirmed by Maplesoft support (probably 2 years ago) and it may be fixed in Maple 2015 (and even Maple 18) but at present none of my machines can run higher than Maple 17 (and no, I will not or cannot upgrade to OSX beyond 10.6.8). I am a bit hesitant to try overwriting the plots:-display package (e.g. using the one from Maple 16 that does not show the bug) until I have at least some understanding what is going on under the hood, since I do not want to kill my Maple installation(s).

Anyone knows how I can inspect plots:-display?

M.D.

@Carl Love @tomleslie

Yes, the DynamicSystems package would require me to recode the whole simulation (which is what produces my data). That may be of advantage eventually, but the whole model I have is suffuciciently complicated that such a recoding would be by no means a trivial task.

But I did investigate the arrayplot in some more detail and it turns out it is useful in that it allows me at least to get the X axes to line up. It is not really what I want it to be, but it gets closer than what I have been using so far, so there is some progress. I played around with plots:-dualaxisplot and so far did not like the results, but there is probably more I can try with it. I also looked at the plot structures documentation but I do not quite see how I can program these to produce a plot closer to my liking. Maybe this is a little beyond what Maple can deliver right now, as far as plotting is concerned...

Thanks much for looking into this,

M.D.

@Carl Love You are of course correct. Don't know what I was thinking...

It would be interesting to know what the OP gets if he just loads CUDA.

On second thought, maybe his libname got messed up? If that gets overridden e.g. by assigning his private library area then Maple likely cannot find the system packages anymore.

M.D.

@Carl Love Ok, I understand that specfunc({sum,Sum}) selects an expression of the form sum() or Sum(). I do not really see why/how patfunc(`+`,anything) applies to the argument of the sum or Sum, i.e. restricts selection to expressions like Sum(A+...+B). I have not actually tried this; but off-hand I would suspect the And never can be True… (implying that Sum() is not of type `+` which I think is the case).

M.D.

@Carl Love Uhh, you got me there... It actually did work in my sheet; but only because I had a local definition of your original function masking my version, which is in a library. Since I really needed the result, first and foremost, I did not double-check this. The problem of course is the factor 2 in front of the sums.

Now, my DistSum still works in my case when I map it over the expression, so all is not lost (and I got where I needed to get already). I'll study your modification to see if I can figure out what exactly it does. I tend to get queasy if I just adopt something without understanding it so I need to do that before it becomes a part of my library. I know the help page you refer to, and yes, it is dense, and I won't claim to understand it.

Many thanks,

M.D.

@Axel Vogt 

@Carl Love

@acer

Thanks much all three of you; all answers/comments were helpful.

The suggestion to use the inert Sum rather than sum is absolutely correct; I had already gotten there (after posting my question) but that in itself did not solve my problem. I actually did check SumTools but did not find anything in there that looked like it would help. It seems like it should have such (rather elementary but useful) functions.

I voted Carl's solution best since it is terse but quite understandable, and I was able to easily modify it to add some checks for validity (as it was; it would merrily distribute the Sum across the factors of a product...No!). This in no way is meant to disregard or disparage Axel's solution. Incidentally, I do not use a range here as the whole computation is done in a general context where the range is not known (and not relevant either). Application of the formulae generated of course will involve ranges.

In case someone is interested, here is the "final" version of Carl's procedure with some added validity checks; coded as a function. Casual testing implies that it works at least for simple cases, and it can be mapped over the terms of simple expressions. It is still a one-liner although I break these up to avoid awkward line breaks.

DistSum:=(xpr)-> `if`(op(0,xpr)='Sum' or op(0,xpr)='sum',\
                                 `if`(type(op(1,xpr),`+`),\
                                       subsindets[flat](xpr, specfunc({sum,Sum}),\
                                                               xpr-> map(op(0,xpr), op(xpr))),\
                                       xpr),\
                                 xpr);

Thanks much,

M.D.

@jbail I just copy-pasted the commands into Maple 15 (fresh sheet, restart on top) and get a plot like Kitonum's.

Are you by any chance in Document mode (i.e. typeset input)? You should be in Worksheet mode whe you paste (red Courier font, Maple prompt at begining of a line.

M.D.

Well, even if it does not print, Maple still handles the expression and you can work with it. Simplify often expands along the way which may cause the long expression. You can try all the usual tricks of the trade to reduce this (I often like collect); and sometimes it helps adding an assumption if you know the domain of some of your variables.

M.D.

@Carl Love Interesting, the use of Record in your example. I had not thought about that.

Below is what I came up with for the exponential averager I need. I want it to return the average after each sample added. I don't need the ability to give a list as argument---I can always map the add function over the list---but I do need a reset function.

The nice thing about doing it this way---besides running several averagers in parallel---is that I can replace the exp. averager with e.g. a douple exp. averager without changing the code that uses it at all, just need to change the constructor. This is one of the features in Maple I cherish.

Mac Dude

--------------

ExpAvg:=proc(lambdap)
 
  return module()
  export add,samples,reset;
  local lambda:=lambdap;
  local avg:=0;
  local n:=0;

    add:=proc(sample) # add sample, return average
      local i;
      avg:=`if`(n>0,lambda*sample+(1-lambda)*avg,sample);
      n:=n+1;
      return avg;
    end proc; # add

    samples:=proc() # return number of samples
      return n
    end proc; # samples

    reset:=proc() # reset this averager
      n:=0;
      avg:=0;
    end proc; # reset

  end module; # the generated module

end proc: # the constructor

ExpAvg.mw

@tomleslie Yeah, I was thinking along those lines initially, too. But I want this to be general & would be concerned tht it is a bit fragile a construct. So I think the module factory is the way to go.

Thanks,

M.D.

@bfathi Ok, without orbit.sav I cannot really run the code. G and Mz are undefined; I assume orbit.sav sets these as well as the arrays XS and YS.

For me the code gets stuck in the loops.

One observation: You use arrays X and Y without declaring them beforehand. While this can be done it leaves Maple to decide what kind of array to use, and I have found I disagreed with Maple's decision. I suggest you do something like
X:=Vector(329,datatype=float);

and the equivalent for any other array you plan to use. With the datatype you also force an evalf and avoid storing symbolic expressions of ever increasing length. If the evalf cannot fully proceed you will get an error and know where the code goes wrong.

Also note that you want to set the starting point for the loop involving n. Not a big deal; but better to know exactly what is being done.

I don't have time right now to continue but hopefully can get back to this later, in the evening.

M.D.

It may not be clear what you mean by "corner". The diagram suggests you are referring to the angles against the vertical axis (what I guess you call the "neutral axis").

I'd first see how to get the angles between the vectors a, b and c from the "corner" angles and the known values of k. That part you can debug by itself. Once you have those, you will be able to solve your problem. It is likely that you get multiple solutions, and you need to settle on the one(s) you want in order to be able to continue.

As far as solving the equations in your worksheet (whatever they mean); you may be able to make progress by replacing the cosine terms with constants, trying to solve, and putting the cosine terms back into the solution. It is probably easier using freeze and thaw for that (rather than subs); or use frontend (which I am having a very hard time with, so I rarely use it). Check the docs for these; it is a bit too detailed to explain here and in any case I would need to fiddle with it to get it right (& I really do not have the time for that right now). It is better you do this on your own.

Yes, there are several obvious errors.

"whit(plots)" on line 1 should be "with(plots)"

Further down there are a number of semicolons (;) missing.

I'd watch out for the
i:='i'; j:=i+1;

construct and make sure it does what you want. If it fails, consider defining j as a function:

j:=ix-> ix+1;

and then reference j(i) when you use it. Note that I deliberately use another name as argument to make clear what this is.

You need to do some basic checking on your code. While Maple's error messages can be cryptic, these ones you should be able to find yourself. Do note that, if Maple returns a function unevaluated (e.g. returns "whit(plots)") it means the function is not (yet) defined; either by design or by error.

M.D.