acer

32358 Reputation

29 Badges

19 years, 332 days
Ontario, Canada

Social Networks and Content at Maplesoft.com

MaplePrimes Activity


These are answers submitted by acer

The 3rd suggestion below, of faking the axes on the smaller, acaled & shifted subplot, is the kind of thing that some people enjoy coding.

The Label in the worksheet's 1st suggestion doesn't inline properly to this site's maplenet, so here that one is (after manually doing bmp-jpg just so I could upload it). I didn't bother inlining the density plot from the 2nd suggestion below, as it makes the worksheet huge.


 

restart:

with(ImageTools):

1st way, by printing both plots to image file, reading in and overlaying those images,

and placing resulting image on a Label component. This is easy.

plotsetup(default):

p1:=plot(sin(x),x=0..2*Pi,tickmarks=[decimalticks,default],
         axesfont=[Helvetica,10,bold],color=blue):

fn:=cat(kernelopts(homedir),"/fn.bmp"):

plotsetup(bmp,plotoutput=fn,plotoptions="height=400,width=600");

print(p1):
plotsetup(default):

fclose(fn): Threads:-Sleep(1):

img:=Read(fn):
rtable_dims(%);

p2:=plot(cos(x),x=0..2*Pi,tickmarks=[decimalticks,default],
         axesfont=[Helvetica,20],labels=[``,``],color=blue):

fn2:=cat(kernelopts(homedir),"/fn2.bmp"):

plotsetup(bmp,plotoutput=fn2,plotoptions="height=400,width=600");

print(p2):
plotsetup(default):

fclose(fn2): Threads:-Sleep(1):

img2:=Read(fn2):

img2:=Scale(img2,0.40):
img2:=Threshold(img2,0.85):
rtable_dims(img2);

#View(img2);
#View(imgnew);

imgnew:=SetSubImage( img, img2, 20, 330 ):
Write(fn,imgnew):

DocumentTools:-SetProperty(Label0,image,fn);

 
 

 2nd way, as density plot produced from overlaid images.

This would look better it the axes were thickened up, and perhaps the curves too, before the earlier

stage of conversion to image. Also GUI performance of 2D densitiy plots is very slow. Don't right-click on it.
nb. This would be also easier if one could place a background image on a regular plot.

 

restart:

with(ImageTools):

plotsetup(default):

p1:=plot(sin(x),x=0..2*Pi,tickmarks=[decimalticks,default],
         axesfont=[Helvetica,16],color=blue,thickness=3,
         axis=[thickness=3]):

fn:=cat(kernelopts(homedir),"/fn.bmp"):

plotsetup(bmp,plotoutput=fn,plotoptions="height=400,width=600");

print(p1): plotsetup(default):

fclose(fn): Threads:-Sleep(1): img:=Read(fn):

p2:=plot(cos(x),x=0..2*Pi,tickmarks=[decimalticks,default],
         axesfont=[Helvetica,30],labels=[``,``],color=blue,
         thickness=5,axis=[thickness=4]):

fn2:=cat(kernelopts(homedir),"/fn2.bmp"):

plotsetup(bmp,plotoutput=fn2,plotoptions="height=400,width=600");

print(p2): plotsetup(default):

fclose(fn2): Threads:-Sleep(1): img2:=Read(fn2):

img2:=Scale(img2,0.40):
img2:=Threshold(img2,0.85):

imgnew:=SetSubImage( img, img2, 20, 330 ):

m,n := ImageTools:-Width(imgnew),ImageTools:-Height(imgnew):

C:=ArrayTools:-Alias(imgnew,[m*n,3]):

ArrayTools:-DataTranspose(C,n,m,1):

P:=ArrayTools:-Alias(C,[m,n,3]):

P:=ArrayTools:-FlipDimension(P,2):

p:=plots:-densityplot(1,x=1..m,y=1..n,axes=none,scaling=constrained,

                      style=patchnogrid,gridlines=false,grid=[m,n]):

p:=subsindets(p,specfunc(anything,COLOR),z->COLOR(RGB,P)):

p;

 

 3rd way, by faking the axes on the subplot as a combination of lines and text, and then using plottools
to transform it -- scale and shift -- and then using plots:-display to overlay in the usual way. Left as

exercise for the reader, because it is finicky work to get the faked axes to look exactly the same as usual.

 

 

Download plotinplot2.mw

acer

I'm glad that people new to Maple don't have to figure out existing code with statements such as `yrt`, `corp`, `esu`, and `eludom`. I don't see that much merit in being terse and slick at the expense of clarity.

These choices were made apparent in Maple 6 (...in Jacques' time), and this all reflects the preferences of the designers and decision makers.

acer

Suppose that you have implemented a Maple program which can be easily split (parellelized) at a high juncture in the code, but for which each computational "thread" involves thread-unsafe operations (eg. calls to Library routines such a `int`, `is`, `limit`, etc). You may not be able to apply the Threads package in order to parellelize the total computation, because the thread-unsafe procedures could clobber each others' global state. In such a situation, you might be able to gain a performance improvement on a single multicore machine by using the Grid toolbox.

Perhaps someone could find the time to illustrate with an example, showing not only a performance gain from using Grid but also demonstrating failure via Theads. The reason I mention parellelizing at a high point in the code is that overhead of Grid is probably not negligible. I'm just guessing, but perhaps a difficult symbolic definite integration example using several methods in parallel.

acer

If you keep having these kinds of troubles (this is the thrid time, at least?~), then why not switch from 2D Math input to 1D Maple notation?

acer

I found it tough to get Maple to look at y=7, x=1..3 as one of the canditate ranges for which the objective is purely real. I was trying to code it in such a way as to give Maple as little help and prodding as I could.

One thing that made this tougher than it ought to be is that I wasn't sure how hard I should try and get Maple to go from something like, {x>=1,x<=3} to x=1..3.  (See the member of newS in the attached worksheet, which correspondes to this.) This aspect makes a difference because, Optimization:-Maximize seems to go awry when this restriction is passed as a pair of constraints, but succeeds when it is passed as simple bounds. See end of attached worksheet.

 

restart:

interface(warnlevel=0):

obj:=sqrt((x-1)*(y-x))+sqrt((7-y)*(1-x))+sqrt((x-y)*(y-7));

((x-1)*(y-x))^(1/2)+((7-y)*(1-x))^(1/2)+((x-y)*(y-7))^(1/2)

plot3d([Re,Im](obj), x=-2..3, y=0..11, color=[red,blue], axes=box);

S:=[solve(evalc(Im(obj)))]:  # takes a little while

newS:={seq(solve(s union {x>=-2,x<=3,y>=0,y<=11}),s in S)};

{{x = 1, y = 1}, {x = 1, y = 4}, {x = 1, y = 7}, {x = 1, 1 <= y, y <= 7}, {x = y, 0 <= y, y <= 1}, {y = 7, x <= 3, 1 < x}}

for s in newS do
   try
      seql,sineq:=selectremove(type,s,`=`);
      ans:=Optimization:-Maximize(eval(obj,seql), sineq,
                                  initialpoint=[x=1/2,y=1/2]);
      print(ans[1],[ans[2][],seql[]]);
   catch:
   end try:
end do;

2.64575131106459072, [y = HFloat(5.551115123125783e-17), x = y]

for s in newS do
   try
      seql,sineq:=selectremove(type,s,`=`);
      ans:=Optimization:-Maximize(eval(obj,seql), sineq,
                                       initialpoint=[x=1,y=5]);
      print(ans[1],[ans[2][],seql[]]);
   catch:
   end try:
end do;

3., [y = HFloat(4.000000000000039), x = 1]

-0., [y = HFloat(1.0), x = y]

#plots:-display(Array([
#   plot([Re,Im](eval(obj,x=1)),y=0..11,view=0..5),
#   plot([Re,Im](eval(obj,y=7)),x=-2..3,view=0..5),  # method above missed this one
#   plot([Re,Im](eval(obj,x=y)),y=0..1,view=0..5) ]));

 

Download toughopt1.mw

For some reason maplenet is not inlining the resst of the sheet. See the link for a bit more.

acer

Is it all T's except at the "square" positions? Let's consider H=0 and T=1, (This could be quite wrong, I've only had one coffee.)

> restart:   
                               
> seq((numtheory[tau](i)-1 mod 2),i=1..50); 

0, 1, 1, 0, 1, 1, 1, 1, 0, 1, 1, 1, 1, 1, 1, 0, 1, 1, 1, 1, 1, 1, 1, 1, 0, 1, 1,

    1, 1, 1, 1, 1, 1, 1, 1, 0, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 0, 1

> evalf[30](sum((1/2^i)*'(numtheory[tau](i)-1 mod 2)',i=1..1000));

                        0.435531586394061420665270725728

> evalf[30](1-sum(1/2^(i^2),i=1..infinity));                      

                        0.435531586394061420665270725728

acer

So, your `eq` is just an equality constraint.

restart;

eq := x = 11/500+6*u*(1/25)-u^2
          +surd(-174200000000*x^5+13045131*E^4*u^4,5)/(65225655*E^4*u^4):

E:=0.1;
ans:=Optimization:-Maximize(x,{eq},u=0..1,x=0..1,initialpoint=[x=1,u=0.5]);
eval((rhs-lhs)(eq),ans[2]);
E:='E':

sols:=seq([E,Optimization:-Maximize(x,{eq},u=0..1,x=0..1,
                                    initialpoint=[x=1,u=0.5])],
          E=0.1..1,0.1);

seq( eval((rhs-lhs)(eq),[s[2,2][],E=s[1]]), s in sols );

acer

Try it with evalf(M) as the first argument to SingularValues.

(This is just an easily fixable bug.)

acer

This appears to be another instance of `sum` versus `add`.

acer

The solutions you've found are all "pretty good" fits to the data. It's possible to obtain a close fit with Optimization:-Minimzie rather than NonlinearFit, but the optimality tolerance has to be tightened.

Also, the solution space has more than just one local minimum, which makes it tough for a local optimizer such as Optimization:-Minimize which will return whichever local minimum it converges to (and not necessarily the globally optimal of all such local minima).

I've created a univariate objective, which fortunately happens to be pretty smooth for your example. But in general there's no reason for it to be so well behaved. (I feel that I got lucky, with the tricks.)

This is just the kind of example where I would expect a curvefitting routine based on global optimization techniques (such as DirectSearch offers) to do much better by default.

 

restart; with(Statistics)

Edata := [-.462124699683824428309502010103, -.527751016544377322226200091423, -.598698303952974890384265763998, -.674809251210729440119324699197, -.75601431564264356773357700158, -.84227676885295360042125209345, -.93357527229747985234963542418, -1.02989666231439422969552894410, -1.13123248073042925355643641564]:

Ez := proc (Z, Zc) options operator, arrow; -(1/2)*Zc^2-a*(Z-Zc)-b*(Z-Zc)^(3/2)-c*(Z-Zc)^2-d*(Z-Zc)^(5/2)+e*(Z-Zc)^3+f*(Z-Zc)^(7/2) end proc:

Ezpw := proc (Z, Zcc) options operator, arrow; piecewise(Zcc < Z, Ez(Z, Zcc)) end proc:

obj := add((eval(Ezpw(Z, Zc), Z = Zdata[i])-Edata[i])^2, i = 1 .. nops(Zdata)):

``

best := infinity:

Optimization:-Minimize(G, .9 .. .95);

[3.59737560152029*10^(-16), Vector(1, {(1) = .910849324507901})]

sol;

[0.359738773129374952e-15, [Zc = HFloat(0.9108493350527227), a = HFloat(1.1425316489552044), b = HFloat(0.1742451497775256), c = HFloat(0.7695947603372874), d = HFloat(0.14059245764891992), e = HFloat(0.023045433394426117), f = HFloat(0.008513627717854673)]]

plots:-display(Array([plot(G, .900 .. .95, view = 0 .. 0.1e-5, style = line), plot(G, .900 .. .93, view = 0 .. 0.1e-7, style = line), plot(G, .900 .. .92, view = 0 .. 0.1e-8, style = line)]))

 

Download maphelp_modif.mw

acer

I would disagree that the natural "factor" in your example is 1/3.

The current posting relates to a recent Question here. Alejandro mentioned that using ``() to block automatic simplification had issues because the reversing operation `expand` could hit too much. But this is also true of prefix `%` and `value`, since applying `value` can meddle with other inert subexpressions such as calls to `Int`, etc. The essential problem stems from trying to make a mechanism serve two purposes. It so often turns out in the long run to be a bad idea whenever anyone tries to hijack a pre-existing facility and try and make it do double-duty. So below I'll avoid both those mentioned ways, and try and use some made-up wrapping procedure (the concept of which has also been done a few times before on this site). As name for the wrapping and unwrapping functions one can try and choose names that aren't already used for other purposes.

restart:

`print/_`:=proc(t) ``(t); end proc:

block:=proc(expr) local c,r;
           c:=content(expr);
           r:=expr/c;
           _(c)*r;
       end proc:

unblock:=proc(expr)
             subsindets(expr,'specfunc(anything,_)',t->op(t));
         end proc:

a:= (1/3)*x^2 + (1/6)*x - (1/12);

                        1  2   1     1 
                        - x  + - x - --
                        3      6     12

A:=block(a);

                      /1 \ /   2          \
                      |--| \4 x  + 2 x - 1/
                      \12/   
              
unblock(A);

                        1  2   1     1 
                        - x  + - x - --
                        3      6     12

acer

Procedures are of type last_name_eval. (Visit that link and see near the end for a tip on returning a procedure as the result from another procedure.)

> Phi := unapply(A1[1], x1, x2);

                           Phi := (x1, x2) -> A1[1]

> Phi;

                                      Phi

> eval(Phi,1);

                               (x1, x2) -> A1[1]

> type(Phi,last_name_eval);

                                     true

acer

Try it with the exponent ``(1/2) instead of 0.5 where `` is two single left-quotes (name-quotes).

If it gets the behaviour you are after, then try applying the expand command to the entire expression and see if that gets you back to normal. (If it works then the sqrt(2) should pop out from the rest.)

acer

If your real, floating-point data is symmetric then try,

B := Matrix(..., shape=symmetric);

Eigenvalues(B);

so as to get an appropriate solver which returns a real result. The argument to `Matrix` above (where I had ... ellipsis) could be an Array, or listlist, or another Matrix the previous version of B, etc.

It's difficult to answer about that error message (indicating a missing argument) without seeing your code. Could you upload a self-contained worksheet that demonstated the problem?

acer

There certainly would be good uses for an "Interpolate" command which, given data, could return an interpolating procedure. Such a returned interpolating procedure might even be made clever in that it could know how to apply the integrated version of its own simple interpolatory formulae (splines, say). But such a thing does not yet exist, as described.

As you say, one need is for very low overhead in on-the-fly interpolation at a given (single) point, which might be done repeatedly but cannot all be done in bulk at once.

Having said that, it might still be possible to speed up your application a bit. If you don't have a great deal of data points then you might be able to use BSplineCurve and instead get an interpolating piecewise polynomial (which could be integrated, once). Or, if sticking with ArrayInterpolation then a little relief might be had with reusable 1-element float[8] Arrays for input/output. Or, possibly, for a fixed degree spline use in ArrayInterpolation then a particular differencing scheme might be able to quickly integrate "exactly".

Any chance you could upload a self-contained worksheet that implemented your current methodology?

acer

First 255 256 257 258 259 260 261 Last Page 257 of 336