## 30305 Reputation

18 years, 256 days

## one way...

It seems that it can be done -- at least for smaller n -- by temporary substitutions a[i]=S[i]+b[i], then resubstituting for S[i] at the end.

The timing increases with n. Perhaps over badly. I lost patience waiting for P(7) to finish. There are also memory consumption issues involved. The OP did mention "small m" in the Question's title, however.

Of course, this is somewhat ad hoc and tailored for the problem at hand. I started by considering replacing a patterned collection of the arguments to the trig calls present. I tried a few, before success.

It might, however, suggest a broader approach. I mean, something other than scanning around speculatively for substitutions or selective freezing.

 > restart;
 > sineExpr:=(m::posint)-> local t,j;                         add(mul(ifelse(j<>t,                                        ':-sin'(a[j]-b[t])/':-sin'(b[j]-b[t]),                                        ':-sin'(a[t]-b[t])),j=1..m),t=1..m):

 > P:=proc(n) local i,e,r,u,S;   e:=sineExpr(n):   r:=[seq(a[i]=S[i]+b[i],i=1..n)];   u:=[seq(S[i]=a[i]-b[i],i=1..n)];   eval(combine(simplify(expand(eval(e,r)))),u); end proc:

 > seq(CodeTools:-Usage(print(P(i))), i=1..6);

memory used=24.02MiB, alloc change=38.01MiB, cpu time=229.00ms, real time=229.00ms, gc time=33.99ms

memory used=4.63MiB, alloc change=0 bytes, cpu time=54.00ms, real time=54.00ms, gc time=0ns

memory used=62.38MiB, alloc change=70.00MiB, cpu time=680.00ms, real time=649.00ms, gc time=135.74ms

memory used=21.23MiB, alloc change=0 bytes, cpu time=192.00ms, real time=193.00ms, gc time=0ns

memory used=124.24MiB, alloc change=-4.00MiB, cpu time=1.25s, real time=1.14s, gc time=230.99ms

memory used=1.58GiB, alloc change=163.77MiB, cpu time=15.36s, real time=13.51s, gc time=3.41s

And further words of caution. Some similar substitutions only work (here, with these particular intermediate steps) only for limited n values, eg.

 > restart;
 > sineExpr:=(m::posint)-> local t,j;                         add(mul(ifelse(j<>t,                                        ':-sin'(a[j]-b[t])/':-sin'(b[j]-b[t]),                                        ':-sin'(a[t]-b[t])),j=1..m),t=1..m):

 > P2:=proc(n) local i,e,r,T;   e:=sineExpr(n):   r:=add(b[i],i=1..n)=T;   eval(combine(simplify(expand(simplify(e,{r})))),(rhs=lhs)(r)); end proc:

 > seq(CodeTools:-Usage(print(P2(i))), i=1..4);

memory used=28.61MiB, alloc change=71.01MiB, cpu time=262.00ms, real time=333.00ms, gc time=0ns

memory used=24.48MiB, alloc change=37.00MiB, cpu time=272.00ms, real time=444.00ms, gc time=80.57ms

memory used=120.70MiB, alloc change=28.89MiB, cpu time=1.26s, real time=1.38s, gc time=245.30ms

memory used=18.29GiB, alloc change=120.27MiB, cpu time=3.66m, real time=3.64m, gc time=55.68s

With an extra step you can obtain both any custom gridline color as well as the nicer curve (directly using the plot command, as you did) in your Maple 2023, as demonstrated in the following attachment.

 >
 >
 >

 >

Unfortunately your color adjustment to gridlines in the axis option does not work directly with plot's default of adaptive=geometric in Maple 2023.

That makes it fall back to adaptive=true, which produces that suboptimal curve portion that you've shown.

That's been fixed in Maple 2024.0, and your example works ok there.

ps. Supplying only gridlines=true produces a slightly darker shade than your forced color=gray.

## Matrix, the column->list...

You could import the whole thing as a Matrix, and then convert the second column to a list.

Eg,

M:=ImportMatrix("example.txt",source=delimited,delimiter=" ");
convert(M[..,2],list);

ps. I suspect that you could also do that first step as,
ImportMatrix("dat_ex.txt",source=Matlab)

And here is some timing comparison. I've used Maple 2023.2.1 here instead, because in my Linux 2024.0 the Iterator's internal fails to Compile (to faster code) and without compilation that Iterator approach takes very much longer.

I've changed the upper ends of the ranges from n down to n-2 in Ronan's approach using piecewise (since i,j,k are each at least 1). That improves the timing to compute all values n=1..30 from 1.872sec down to 1.422sec.

I've also included a further revision to that (end-point) modification of Ronan's approach, in which `if` is used instead of piecewise. Using piecewise here is a poor idea, in terms of performance.

 > kernelopts(version);

Warning, (in s) `l` is implicitly declared local

 > restart; str:=time[real](): s := 0: for n to 30 do     comps := Iterator:-Composition(n, parts = 3);     s := s + add(1/mul(comp), comp in comps);     #print(s); end do: time[real]()-str;

## site problems...

This site has been having several (new) issues in the past month or so.

1) It is slow to access, and ofter very slow
2) sometimes submissions are very slow (and if previewed, sometimes so slow that it seems stuck)
3) the response counter is much more broken than before

Several people have made posts about all this. I don't know whether anyone's fixing it.

## simplify@factor...

 > expr:=(10*(5+sqrt(41)))/(sqrt(70+10*sqrt(41))*sqrt(130+10*sqrt(41))):

 > simplify(factor(expr))

 I see now that Thomas had already mentioned evala@factor. which is great. Also possible is radnormal@factor here.

The following variants involve three calls, without using factor.

simplify(evala( combine(expr) ));

That leads to noticing that the following expression (obtained from combine(expr) in M2024.0 at least) does not produce the same result (or normal form) under repeated calls to either evala or radnormal. This deserves a separate bug report, which I shall make.

 > restart;
 > foo := (5+41^(1/2))*(1/(7+41^(1/2))/(13+41^(1/2)))^(1/2);

 > evala(foo);

 > evala(evala(foo));

## indexing...

I wouldn't prefer the LinearAlgebra:-Column command for any example for which Matrix-indexing is convenient and terse. That's just a personal preference.

Eg,

 > M := LinearAlgebra:-RandomMatrix(5);

 > seq(M[..,i],i=1..5);

## syntax mistakes...

Your code has several invalid calls to int, and several missing commas in the code of the plot calls.

ode_Plots_error_ac.mw

## example...

Do you mean something like the following? We can get similar effects using a loop, or using seq.

 > M := Matrix(4,4): for i from 1 to 4 do   for j from 1 to i do     M[i,j] := [i,j];   end do; end do: M;

 > seq( seq( [i,j], j=1..i ), i=1..4);

## various changes...

Ok, so there are a few issues to deal with.

1) Instead of using _rest to pass along any extra options to the plotting commands I set it up to use a new keyword option plopts. That avoids a clash with using the lowercase names for Qdim's keyword options. There are several stock Library commands that use such an approach. (You can back that out, and use the capitalized names or other, if you prefer. I just wanted to show this to you.)
2) I also made the point and line options work, by having Qdim explicitly name a terse substitute pt for plottools in its uses syntax. I also did the same for plots in that uses. This allows for reasonably terse (but understandable) code.
3) I also made Qdim's lowercase print option work, by invoking the global command by :-print instead.
4) I also put uneval quotes around keyword options in use in Qdim, and uneval quotes around the global versions of (stock unprotected) names used as option values. This is to guard against Qdim breaking if any of these names got assigned some values at the top-level. I used maplemint to discover such.

I suggest you stress test (as your examples do).

 > restart;
 > # changed  points   to point #          linetype to line #          scl      to scale   # lower case dont know if clashes #          clr      to colour #          prnt     to print
 > Prntmsg:="y"; Geomclr:="Blue";

 > Qdim := proc(P1, P2, {Q:=[ NULL,:-align={':-left'}]},                      {vec1::list:=[NULL,:-align={':-below',':-right'}]},                      {scale::{list,numeric}:=1},                      {leader::{-1,0,1}:= 1},                      {dimofset::{1,2,3,4,5,6}:=1},                      {point::list:=['color' = ':-blue', symbol = ':-solidcircle', 'symbolsize' = 8]},                      {line::list:=['thickness'=3]},                      {colour::`string`:="blonde"},                      {print:=Prntmsg},                      {plopts::list:=[]}) description " plots Quadrance symbol and value (Q)"; uses pl=plots, pt=plottools; # RationalTrigonometry; local a,f,g,h,v1,delta,mp, BoxQ,Qsymbol,pt1,pt2,pt3,pt4,Scl,       ptxtp1,ptxtp2,txtplt1,txtplt2,plttyp,thk,l12,pts12; if P1[1]::{list, Vector[row]} then    pt1:=`if`(P1[1]::Vector[row],convert(P1[1],list),P1[1]);    plttyp:=`if`(nops(pt1)=3,3,2); # print(pt1,plttyp);    ptxtp1:=[op(pt1),op(P1[2..-1])];    txtplt1:=true;    else    pt1:=`if`(P1::Vector[row],convert(P1,list),P1);#print(pt1);    plttyp:=nops(pt1);    txtplt1:=false; end if;     if P2[1]::{list, Vector[row]} then    pt2:=`if`(P2[1]::Vector[row],convert(P2[1],list),P2[1]);#print(pt2);    if plttyp<>nops(pt2) then error `different dimension in inputs`end if;    ptxtp2:=[op(pt2),op(P2[2..-1])];    txtplt2:=true;  else    pt2:=`if`(P2::Vector[row],convert(P2,list),P2);#print(pt2);    if plttyp<>nops(pt2) then error `different dimension in inputs`end if;    txtplt2:=false; end if; #if scl::numeric then     Scl:=scale;  #else   # pt3:=`if`(scl[1]::Vector[row],convert(scl[1],list),scl[1]);    #pt4:=`if`(scl[2]::Vector[row],convert(scl[2],list),scl[2]);     #if plttyp = 2 then      #  Scl:=evalf(sqrt(Quadrance(pt3,pt4,"b","n")/Quadrance(pt1,pt2,"b","n")));      #else       # Scl:=evalf(sqrt(Quadrance(pt3,pt4,"n")/Quadrance(pt1,pt2,"n")));    # end if; #end if; if plttyp =2 then      a :=Scl* dimofset*0.1*leader*[-pt2[2]+pt1[2], pt2[1]-pt1[1]]; #dimension leader lines   else      a :=Scl* dimofset*0.1*leader*[-pt2[2]+pt1[2], pt2[1]-pt1[1],pt2[3]-pt1[3]]; #dimension leader lines end if; delta:=Scl*(pt2-pt1); mp:=(1-1/2)*pt1+1/2*pt2;  #midpoint of line #print(Scl,a,delta,mp); if plttyp = 2 then     BoxQ:=[[a[1]+mp[1]+.05*delta[1]-.01*delta[2],a[2]+mp[2]+.05*delta[2]+.01*delta[1]],            [a[1]+mp[1]-.05*delta[1]-.01*delta[2],a[2]+mp[2]-.05*delta[2]+.01*delta[1]],            [a[1]+mp[1]-.05*delta[1]+.01*delta[2],a[2]+mp[2]-.05*delta[2]-.01*delta[1]],            [a[1]+mp[1]+.05*delta[1]+.01*delta[2],a[2]+mp[2]+.05*delta[2]-.01*delta[1]]];  else     BoxQ:=[[a[1]+mp[1]+.05*delta[1]-.01*delta[2],a[2]+mp[2]+.05*delta[2]+.01*delta[1],a[3]+mp[3]+.05*delta[3]],            [a[1]+mp[1]-.05*delta[1]-.01*delta[2],a[2]+mp[2]-.05*delta[2]+.01*delta[1],a[3]+mp[3]-.05*delta[3]],            [a[1]+mp[1]-.05*delta[1]+.01*delta[2],a[2]+mp[2]-.05*delta[2]-.01*delta[1],a[3]+mp[3]-.05*delta[3]],            [a[1]+mp[1]+.05*delta[1]+.01*delta[2],a[2]+mp[2]+.05*delta[2]-.01*delta[1],a[3]+mp[3]+.05*delta[3]]]; end if; if leader<>0  then    thk:='thickness'=0  elif line<>[] then    thk:=line[op(ListTools:-Search(op(select(has,line,'thickness')),line))];  else    thk:='thickness'=0; end if; Qsymbol:= [pt:-line(pt1,pt1+a,'thickness'=0,plopts[]),            pt:-line(pt1+a,mp-.05*delta+a,thk,plopts[]),            pt:-line(op([BoxQ[1],BoxQ[2]]),'thickness'=0,plopts[]),            pt:-line(op([BoxQ[2],BoxQ[3]]),'thickness'=0,plopts[]),            pt:-line(op([BoxQ[3],BoxQ[4]]),'thickness'=0,plopts[]),            pt:-line(op([BoxQ[4],BoxQ[1]]),'thickness'=0,plopts[]),            pt:-line(pt2+a,mp+.05*delta+a,thk,plopts[]),            pt:-line(pt2,pt2+a,'thickness'=0,plopts[])];      if point<>[]then            pts12:=pt:-point([pt1,pt2],point[]);  #end points else pts12:=NULL; end if; if line<>[] then l12:=pt:-line(pt1,pt2,line[]);   else   l12:=NULL;   end if; if plttyp = 2 then     f := pl:-textplot([op((pt1+pt2)/2+1*(a)),op(Q)]);  #Quadrance text    if txtplt1 then       g := pl:-textplot( ptxtp1);     else       g:=NULL;    end if;    if txtplt2 then       h := pl:-textplot( ptxtp2);     else      h:= NULL;    end if;     v1:=pl:-textplot([op((pt1+pt2)/2), op(vec1)]);#`#mover(mi("v"),mo("⇀"))`[1,2]   {':-below',':-right'}  else    f := pl:-textplot3d([op((pt1+pt2)/2+1*(a)),op(Q)]);    if txtplt1 then          g := pl:-textplot3d( ptxtp1);       else          g:=NULL;    end if;    if txtplt2 then          h := pl:-textplot3d( ptxtp2);       else          h:=NULL;    end if;     v1:=pl:-textplot3d([op(pt1+pt2)/2,op(vec1)]);#`#mover(mi("v"),mo("⇀"))`[1,2]   {':-below',':-right'} end if; #print("pltyp" , plttyp)  ; if print="y" then :-print(cat("Quadrance symbols  ",colour,"  geometry")); end if; if leader = 0 then         pl:-display(pts12,f,g,h,Qsymbol,'axes'=':-none','scaling'=':-constrained');  else                 pl:-display(pts12,l12,f,g,h,Qsymbol,'axes'=':-none','scaling'=':-constrained');      end if; end proc:
 > maplemint(Qdim);

Procedure Qdim( P1, P2, { Q := [NULL, :-align = {':-left'}], colour::string :=
"blonde", dimofset::{1, 2, 3, 4, 5, 6} := 1, leader::{-1, 0, 1} := 1,
line::list := [('thickness') = 3], plopts::list := [], point::list :=
[('color') = (':-blue'), symbol = (':-solidcircle'), ('
symbolsize') = 8],
print := Prntmsg, scale::{list, numeric} := 1, vec1::list := [NULL,
:-align = {':-below', ':-right'}] } )
These local variables were never used:  pt3, pt4
These local variables were assigned a value, but otherwise unused:  v1

 > P1:=[2,3,6]: P2:=[1,5,7]: P3:=P1+1/3*(P2-P1)

 > plt1:=Qdim([P1,typeset("P1=",P1),align=below],[P2,P2,align=right ],            dimofset=2,Q=["Q12\n",align =above],plopts=[colour=black]);

 > plt2:=Qdim([P1,typeset("P1=",P1),align=below],[P3,P3,align=right ],            scale=2,Q=["Q13\n",align =above],print="n"); #scl affects the linine up of dimesion lines prnt="n" no message displayed

 > plt3:=Qdim(P3,P2,scale=1,Q=["Q23\n",align =above],point=[],            line=[],colour="RED"); #doesn't plot points and line. shows use of `clr`

 > plots:-display(plt1,plt2,plt3)

 > P5:=[1,2]:P6:=[-1,4]:
 > Qdim([P5,"P5   ",align={below,left}],[P6,"   P6",align={right}],      Q=[typeset("\nQ[5,6]"),align={below,left}],line=[linestyle=dash],      point=[symbolsize=20,symbol=solidcircle,colour=orange]);

 > Qdim([P5,"P5   ",align={below,left}],[P6,"   P6",align={right}],      leader=0,Q=[typeset("\nQ[5,6]"),align={below,left}],      point=[symbolsize=18]);

 >

## one way...

I've used Maple 2024.0 here.

Since simplify is not by itself (well, no longer, since M2023) finding the trig simplification then the key to my approach above is to force the trig combine of the sin and cos terms in the second expression. The subsequent trig expand gets the simpler thing with just cos.

I shall submit bug reports against the regression & weaknesses in simplify.

## plottools:-line, and memory size, etc...

If you're going to use plot3d to construct 3d lines then you might as well pass the option grid=[2,2] and cut down the wasted memory. The default is [49,49], but you only need the end-points, so without that option your plot3d construction is unnecessarily about 25 times too big.

Similarly, if you're using spacecurve for mere lines and want solid colors then you could pass numpoints=2 to that command, again to avoid unnecessarily wasting time and space. (It's hard to measure the time gain, since the example is already very quick.)

That might all seem pedantic, but those kind of considerations can really come into play when one starts animating/exploring plots, etc.

But why not just use plottools:-line, with the two Vectors you already have? I mean, instantiating them both at the end-points for lambda. That get's that memory efficiency directly.

 > restart
 > with(plottools): with(plots,display):
 > l:=([2,-3,1],<3,7/9,6>):   # 3d line point + vector
 > P:=[7,-8,9]:
 > pl:=`+`~(lambda*l[2],l[1]): #3d line as vector eqn
 > vnl:=`-`~(pl,P) : #vector from Point P to 3D line
 > vnl.l[2] assuming `real` ; #dot product of vectors= 0 when perpendicular

 > sol:=solve( {  }, [lambda] )[];

 > intP:=eval(pl,sol):  #intersection point
 > l2:=P,eval(vnl,sol) :  #perpendicular 3D line through P
 > pl2:=`+`~(lambda*l2[2],l2[1]): #3D line as vector eqn
 > display(       line(eval(pl,lambda=-.5),eval(pl,lambda=1.8),thickness=0,colour=orange),       line(eval(pl2,lambda=-.5),eval(pl2,lambda=1.8),thickness=0,colour=purple),       point(P,colour=blue ,symbolsize=15,symbol=solidsphere),       point(l[1],colour=green ,symbolsize=15,symbol=solidsphere),       point(eval(pl,sol),colour=red ,symbolsize=15,symbol=solidsphere),       arrow(l,0.2, 0.4, 0.1,colour=green),       arrow(l2,0.2, 0.4, 0.1,colour=blue),       axes=normal,scaling=constrained);

## possible...

You might be able to get by with,

is(simplify(tmp - r) = 0)

quite often.

## some ideas...

I'm not sure what you want, as "simpler".

I haven't included anything like alias or LargeExpression:-Veil.

But maybe there are some ideas here for you (including a variety of "expression size" metrics),
expression_to_simplify_ideas.mw

## global...

Your call to the parse command constructs the global name `a`, which is not your procedure's assigned local a.

 > restart;