tomleslie

13876 Reputation

20 Badges

15 years, 185 days

MaplePrimes Activity


These are replies submitted by tomleslie

@rsweet 

Use the big green up-arrow in the Mapleprimes toolbar to upload your worksheet.

"Pictures" of code are generally pretty useless.

 

@kfli 

I'm not generally a big fan of writing code in one language and translating it into another. In my experience this (usually) works about 90% of the time. You may be able to get this up to 99% by writing the original code in a particular fashion - whihc usually means "Don't do anything clever".

In your particular case, I can write a procedure which accepts an arbitrary length vector, processes the elements of this vector and returns an indexable quantity (in MapleSpeak, a table) which Matlab will interpret as a vector.

For any Maple purist reading this, I wouldn't do it this way unless I had to translate/export the relevant code to Matlab, so don't give me a hard time about the 'sloppy' Maple code.

Obviously there has to be some systematic/programmatic method of generating elements of the output vector from elements of the input vector (whihc only you will know!).

In the following "toy example, if the input vector is [x1, x2,..xN], then the output vector is [x1+x1^2, x1+x2^2,... x1+xN^2], although any other (programmable combination would be possible). The resulting code *seems* to translate successfully to Matlab

restart;
myV:=proc(x)
                    local N:=numelems(x),
                              j, X;
                    for j from 1 by 1 to N do
                         X[j]:=x[1] +x[j]^2;
                    od:
                    return X;
           end proc;
x:=Vector([x1, x2, x3, x4, x5, x6]):
p:=eval(myV(x));
CodeGeneration:-Matlab( myV);
CodeGeneration:-Matlab( myV,
                                            output="C:/Users/TomLeslie/Desktop/test.m"
                                          );

 

@litun 

If you get problems with "not numbers" - then use the ".csv" option given in my original post.

I can't actually check this in Maple 17 (don't have anything that old) but I have just verified it in

Maple 18
Maple 2015
Maple 2016
Maple 2017
Maple 2018

and it works in all of these with no problems/warnings/issues at all

@Pierre-alexandre 

is to appreciate that you have a lot of data.

Before doing anything else, I would think long and hard about how to organise (ie store) this data, so that I could always easily access whatever data I needed. This might nvolve writing one or more "access" functions.

From reading your descriptions (and I might be wrong!), you seem to have

  1. Different timesteps (snapshots, frames, whatever)
  2. Different types of molecule
  3. Distinct molecules within each moleculeType
  4. For each individual molecule, some "sites" (given by coordinates?), and maybe some other properties - eg some measure of the "size" of each molecule

I would argue that you need to organize data storage+access function(s), so that issuing a command as simple as

access( timestep, moleculeType, moleculeNumber, moleculeProperty)

will always return the desired information for further processing. In the access() function you also have to consider whether to handle ranges for any of the passed parameters, as well as the return type (particularly if you are going to use ranges

##########################################################################

In the above I have implied an access hierarchy,

timestep->moleculeType->moleculeNumber->moleculeProperty

but depending on the caclulations you intend to perform, a different hierarchy may be more appropriate: for example

moleculeType->moleculeNumber->timestep->moleculeProperty

would work equally well. Only you know which hierarchy is more appropriate for the calculations you are planning to perform

##########################################################################

This might sound like a lot of work before you get to do any interesting calculations. Believe me, if you are going to do anything other than a single (fairly trivial) calculation, organizing your data+access function(s) will save a lot of pain in the future!

@acer 

I've seen this type of command (eg 'Typesetting:-mn') presented on this site before (possibly by you!) so I tried to figure out what was being done. The first problem I had, was that Maple help does not list eg 'mn()' as a command in the Typesetting() package. There is a statement in the help, (my emphasis)

Most other Typesetting exports correspond to either internal-use routines or typesetting tags, which roughly correspond to MathML tags. The only one of these that may be useful is the Typeset export,

I happen to know that both 'mn()' and 'mrow()' are used as MathML tags.

exports(Typesetting) shows that these two commands (along with many undocumented others) are available(?) for use. Would it be generally safe to assume that any command returned by exports(Typesetting) which "shares a name" with a MathML tag will do pretty much what the latter would do - or would it be safer to regard it as "suck it and see"?

Do you happen to know if the Maple help for the Typesetting() package is going to be enhanced/extended to cover some/all of these undocumented commands, or are they destined to remain effectively "hidden"

@Pierre-alexandre 

See the attached, which splits the input matrix (whenever the '1728' stuff happens) into a list of matrices. The subsequent list  of matrices is then processed to retian only those terms whihc have C1 (could be anything else you want) in the first column

  restart;

  fileToImport := "C:/Users/TomLeslie/Desktop/test.txt":
  M1:=ImportMatrix( fileToImport):
#
# Get the row numbers where 1728 appears
#
  L:=[seq(`if`(M1[j,1]=1728, j, NULL), j=1..op([1,1], M1))];
#
# Split the input matrix into a list of matrices, bounded
# by where 1728 appears
#
  LM:= [M1[1..L[1]-1, 1], seq(M1[L[j]+1..L[j+1]-1,1], j=1..numelems(L)-1), M1[L[-1]+1..-1,1]];
#
# Strip the "generated by VMD" stuff
#
  F:=vec-> <seq(`if`(vec[j]="generated by VMD", NULL, <vec[j]>), j=1..op(1, vec))>:
  M2:=F~(LM[2..-1]);
#
# Parse each matrix in the list
#
  G:=mat->Matrix(sscanf~(convert(mat, list), "%s%f%f%f")):
  M3:=G~(M2);
#
# Retain only entries matching iden
#
  iden:="C1":
  H1:= mat-><seq(`if`(mat[j,1]=iden, mat[j,..], NULL), j=1..op([1,1],mat))>:
  H2:= mat-><seq(`if`(mat[j,1]=iden, mat[j,2..], NULL), j=1..op([1,1],mat))>:
  H1~(M3);
  H2~(M3);

L := [1, 8, 13]

 

"[[[[]]],[[["generated by VMD"],["C1        -6.932300       12.540000      -20.260000"],["C1        -7.005000       13.914000      -18.431999"],["C2        -8.143600       11.942000      -19.424000"],["C2        -8.236300       12.972000      -18.327000"],["C3        -5.110600       12.270000      -18.719000"]]],[[["generated by VMD"],["C4        -5.390700       14.318000      -20.143000"],["C5        -4.419600       15.242000      -19.368999"],["C6        -2.534000       15.407000      -20.941000"]]],[[["generated by VMD"],["C1        10.458000      -20.974001       -1.023700"],["C1         9.977000      -19.162001        0.231220"]]]]"

 

[Vector(5, {(1) = "C1        -6.932300       12.540000      -20.260000", (2) = "C1        -7.005000       13.914000      -18.431999", (3) = "C2        -8.143600       11.942000      -19.424000", (4) = "C2        -8.236300       12.972000      -18.327000", (5) = "C3        -5.110600       12.270000      -18.719000"}), Vector(3, {(1) = "C4        -5.390700       14.318000      -20.143000", (2) = "C5        -4.419600       15.242000      -19.368999", (3) = "C6        -2.534000       15.407000      -20.941000"}), Vector(2, {(1) = "C1        10.458000      -20.974001       -1.023700", (2) = "C1         9.977000      -19.162001        0.231220"})]

 

[Matrix(5, 4, {(1, 1) = "C1", (1, 2) = -6.932300, (1, 3) = 12.540000, (1, 4) = -20.260000, (2, 1) = "C1", (2, 2) = -7.005000, (2, 3) = 13.914000, (2, 4) = -18.431999, (3, 1) = "C2", (3, 2) = -8.143600, (3, 3) = 11.942000, (3, 4) = -19.424000, (4, 1) = "C2", (4, 2) = -8.236300, (4, 3) = 12.972000, (4, 4) = -18.327000, (5, 1) = "C3", (5, 2) = -5.110600, (5, 3) = 12.270000, (5, 4) = -18.719000}), Matrix(3, 4, {(1, 1) = "C4", (1, 2) = -5.390700, (1, 3) = 14.318000, (1, 4) = -20.143000, (2, 1) = "C5", (2, 2) = -4.419600, (2, 3) = 15.242000, (2, 4) = -19.368999, (3, 1) = "C6", (3, 2) = -2.534000, (3, 3) = 15.407000, (3, 4) = -20.941000}), Matrix(2, 4, {(1, 1) = "C1", (1, 2) = 10.458000, (1, 3) = -20.974001, (1, 4) = -1.023700, (2, 1) = "C1", (2, 2) = 9.977000, (2, 3) = -19.162001, (2, 4) = .231220})]

 

"[[[["C1",-6.932300,12.540000,-20.260000],["C1",-7.005000,13.914000,-18.431999]]],[[[]]],[[["C1",10.458000,-20.974001,-1.023700],["C1",9.977000,-19.162001,0.231220]]]]"

 

[Matrix(%id = 18446744074377452830), Vector[column](%id = 18446744074377452950), Matrix(%id = 18446744074377453310)]

(1)

 

Download 2impData.mw

@Pierre-alexandre 

but it is trivial to count the number of rows in the matrix at varios points - see the following

  restart;
  fileToImport := "C:/Users/TomLeslie/Desktop/test.txt":
  M1:=ImportMatrix( fileToImport);
#
# Number of rows in the input matrix
#
  nR1:=op([1,1],M1);
  M2:=<seq(`if`(M1[j,1]=1728 or M1[j,1]="generated by VMD", NULL, <M1[j,1]>), j=1..op([1,1], M1))>;
#
# Number of rows after clearing out the
# '1728' and "generated by VMD" stuff
#
  nR2:=op([1,1], M2);
  M3:=Matrix(sscanf~(convert(M2, list), "%s%f%f%f"));
  iden:="C1":
  M4:=<seq(`if`(M3[j,1]=iden, M3[j,..], NULL), j=1..op([1,1],M3))>;
  M5:=<seq(`if`(M3[j,1]=iden, M3[j,2..], NULL), j=1..op([1,1],M3))>;
#
# Number of rows in the output matrix
#
  nR3:=op([1,1], M4);

 

@Pierre-alexandre 

I might be tempted to use a decent editor to fix up the data file before importing into Maple.

Alternatively, if I insert the "offending"  information into the test textfile a few times to produce

 1728
generated by VMD
C1        -6.932300       12.540000      -20.260000
C1        -7.005000       13.914000      -18.431999
C2        -8.143600       11.942000      -19.424000
C2        -8.236300       12.972000      -18.327000
C3        -5.110600       12.270000      -18.719000
1728
generated by VMD
C4        -5.390700       14.318000      -20.143000
C5        -4.419600       15.242000      -19.368999
C6        -2.534000       15.407000      -20.941000
1728
generated by VMD
C1        10.458000      -20.974001       -1.023700
C1         9.977000      -19.162001        0.231220

Then the following code will produce the same results as before

restart;
fileToImport := "C:/Users/TomLeslie/Desktop/test.txt":
M1:=ImportMatrix( fileToImport);
M2:=<seq(`if`(M1[j,1]=1728 or M1[j,1]="generated by VMD", NULL, <M1[j,1]>), j=1..op([1,1], M1))>;
M3:=Matrix(sscanf~(convert(M2, list), "%s%f%f%f"));
iden:="C1":
M4:=<seq(`if`(M3[j,1]=iden, M3[j,..], NULL), j=1..op([1,1],M3))>;
M5:=<seq(`if`(M3[j,1]=iden, M3[j,2..], NULL), j=1..op([1,1],M3))>;

 

@Kitonum 

because you have shown how to do it without using the display() command.

To be honest, I would have *expected* the OP's original code to work - list of plots+list of colors: this will work in most(?) plotting commands.

The point I was trying to make was that if a particular plotting command does not accept a partciular plotting option, then that option can often be applied by using the display() command. As the help for the display() command does state (my emphasis)

Additional options as described in the plot/option and plot3d/option help pages may be provided if they are applicable. If the option is one that applies to the entire plot, such as the title or axes style, then the new option overrides any previously specified in the plots contained in P. However, options such as color or linestyle that apply to individual elements like curves or polygons do not affect previously set colors or linestyles.

 

 

@vv 

The only assumptions I made were

  1. all structures are composed of triangles.
  2. each rhombus is two identical triangles, stacked back-to-back
  3. all triangles have the same altitude

It wasn't difficult to adjust my original to fulfil your second point, ie to ensure that the outermost rhombus is in fact a square. It just means calculating triangle altitude a slightly different way. See mosaic2.mw

Your first point, that both sets of quadrilaterals be rhombuses, caused me a 'Necker Cube' moment. Somehow I hadn't really *looked* at the 'unfilled' shapes - I was focussed on the 'filled' ones. I found this trickier to enforce, since one now needs five different triangle altitudes. It is possible to do this whilst still using the geometry() package and cartesian coordinates, but it got *uglier* than I would like, see  mosaic3.mw .

It did convince me that a better approach would probably involve ditching the geometry() package and using polar coordinates!

Code you supply runs with no problems in Maple 2018 - see the attached

restart;
kernelopts(version);
l := (n+m+sum(a[i], i = 1 .. n)+sum(b[j], j = 1 .. m))*ln(alpha)+n*ln(lambda[1])+m*ln(lambda[2])+lambda[1]*(sum(x[i], i = 1 .. n))+lambda[2]*(sum(y[j], j = 1 .. m))-(sum((2+a[i])*ln(exp(lambda[1]*x[i])-1+alpha), i = 1 .. n))-(sum((2+b[j])*ln(exp(lambda[2]*y[j])-1+alpha), j = 1 .. m));

`Maple 2018.1, X86 64 WINDOWS, Jun 8 2018, Build ID 1321769`

 

(n+m+sum(a[i], i = 1 .. n)+sum(b[j], j = 1 .. m))*ln(alpha)+n*ln(lambda[1])+m*ln(lambda[2])+lambda[1]*(sum(x[i], i = 1 .. n))+lambda[2]*(sum(y[j], j = 1 .. m))-(sum((2+a[i])*ln(exp(lambda[1]*x[i])-1+alpha), i = 1 .. n))-(sum((2+b[j])*ln(exp(lambda[2]*y[j])-1+alpha), j = 1 .. m))

(1)

 

Download badsums.mw

I decided to try coding myself - not because there is anything wrong with yours, but

  1. I thought it was the only way (for me) to really appreciate the symmetries of the pattern
  2. I wondered how difficult it might be to do in a completely different way
  3. Also I was curious as to whether it could be done using only "triangle geometry", with no "trigonometry" at all: I failed on this last point, because I always ended up needing 'tan(Pi/m)', where 'm' is the number of pattern repeats in the circle (ie 32 in the real case)

Anyway, FWIW, see the attached

#
# Initialise and set u some parameters
#
  restart;
  with(geometry):
  ri:=10:       # inner radius of tiling pattern
  ro:=40:       # outer radius of tiling pattern
  h:=(ro-ri)/8: # cos there are 8 triangles along a radius
  m:=32:        # cos there are 32 radial patterns.
                # although any other number will work
#
# Define a line along the x-axis, needed for
# reflection purposes
#
  line( l1,
        [ point(lp1, [0,  0]),
          point(lp2, [ro, 0])
        ]
      ):
#
# Compute triangles along a radius. Two superfluous
# triangles will be created, T0 and T9, just because
# it's easier that way, but these will be ignored
# for plotting purposes
#
  for j from 0 by 2 to 8 do
      triangle( cat(T, j ),
                [ point( cat( T, j, A ),
                         [ ri+(j-1)*h,
                           0
                         ]
                       ),
                  point( cat( T, j, B ),
                         [ ri+j*h,
                           -(ri+j*h)*tan(Pi/m)
                         ]
                       ),
                  reflection( cat( T, j, C ),
                              cat( T, j, B ),
                              l1
                            )
                ]
              );
      reflection( cat(T, j+1),
                  cat(T, j),
                  line( l2,
                        [ cat( T, j, C ),
                          cat( T, j, B )
                        ]
                      )          
                );
  od:
#
# Generate the plot of triangles along a single radius
#
  p1:=draw( [ circle( c,
                      [ point(C,0,0),
                        ri-4
                      ]
                    )(filled=true),
              seq( cat(T,j)(filled=true),
                   j=1..8
                 )
            ]
          ):
#
# Rotate the above plot 'm' times, by the increment 2*Pi/m
#
  plots:-display
         ( [ seq
             ( plottools:-rotate( p1,
                                  2*j*Pi/m,
                                  [0, 0]
                                ),
               j=0..m-1
             )
           ],
           size=[800,800]
         );

 

 

Download mosaic.mw

although I have to admit I'm more impressed by the guy who did the original tiling ;-)

Now about my bathroom remodelling.........

@Ham 

The attached shows how to "reverse" the order of the plots and change the vertical spacing between them.

Can't do much about the Tabulate/title issue becuase as far as I can tell, the tabulate command does not accept a 'Title' option. I can only assume that this falls into the "partially-written-software" category: ie the option is accepted, but has not (yet?) been fully implemented, so setting teh option actually does nothing (but doesn't "error")

  restart;
  with(Logic):
  separation:= 1.0:
  TT:=TruthTable( x[1] &xor x[2] &xor x[3] &or x[4], form=MOD2):
  yAxisLabels:= [ColumnLabels(TT)[1..-2][],Q]:
  rows:= upperbound( TT, 1 ):
  columns:= upperbound( TT, 2):
  xAxis:= time:
  lineColours:= [seq(blue,i=1..columns-1), red]:
#
# Original setup
#
  plotOptions:= xAxis=0..rows-1, title = "Timing diagram",titlefont=[times,bold,18],
                color=lineColours,
                tickmarks= [ default,[seq( [j*(1+separation)=0, 1+j*(1+separation)=1][], j=0..columns-1)]],
                axis[1]= [gridlines= [[seq(j=j,j=0..rows)],colour= gray,linestyle= dash]]:
  toPW:=M->piecewise(seq( [xAxis<j, M[j]][], j=1..rows,1)):
  timeValues:= seq( toPW(TT[j])+ (1+separation)*(j-1), j=1..columns,1):
  timePlots:= plot( [timeValues], plotOptions):
  labelValues:= seq( [-1,0.5+j*(1+separation), yAxisLabels[j+1]], j=0..columns-1,1):
  labelPlots:= plots:-textplot([labelValues], font=[times,bold,12]):
  plots:-display(timePlots, labelPlots);
#
# Plots in "reverse" order
#
  timePlots:= plot( ListTools:-Reverse([timeValues]), plotOptions):
  labelPlots:= plots:-textplot(ListTools:-Reverse([labelValues]), font=[times,bold,12]):
  plots:-display(timePlots, labelPlots);
#
# Plots in reverse order with different vertical
# "separations"
#
  separation:=0.5:
  plotOptions:= xAxis=0..rows-1, title = "Timing diagram",titlefont=[times,bold,18],
                color=lineColours,
                tickmarks= [ default,[seq( [j*(1+separation)=0, 1+j*(1+separation)=1][], j=0..columns-1)]],
                axis[1]= [gridlines= [[seq(j=j,j=0..rows)],colour= gray,linestyle= dash]]:
  timeValues:= seq( toPW(TT[j])+ (1+separation)*(j-1), j=1..columns,1):
  timePlots:= plot( ListTools:-Reverse([timeValues]), plotOptions):
  labelValues:= seq( [-1,0.5+j*(1+separation), yAxisLabels[j+1]], j=0..columns-1,1):
  labelPlots:= plots:-textplot(ListTools:-Reverse([labelValues]), font=[times,bold,12]):
  plots:-display(timePlots, labelPlots);

  separation:=0.25:
  plotOptions:= xAxis=0..rows-1, title = "Timing diagram",titlefont=[times,bold,18],
                color=lineColours,
                tickmarks= [ default,[seq( [j*(1+separation)=0, 1+j*(1+separation)=1][], j=0..columns-1)]],
                axis[1]= [gridlines= [[seq(j=j,j=0..rows)],colour= gray,linestyle= dash]]:
  timeValues:= seq( toPW(TT[j])+ (1+separation)*(j-1), j=1..columns,1):
  timePlots:= plot( ListTools:-Reverse([timeValues]), plotOptions):
  labelValues:= seq( [-1,0.5+j*(1+separation), yAxisLabels[j+1]], j=0..columns-1,1):
  labelPlots:= plots:-textplot(ListTools:-Reverse([labelValues]), font=[times,bold,12]):
  plots:-display(timePlots, labelPlots);

  separation:=0.1:
  plotOptions:= xAxis=0..rows-1, title = "Timing diagram",titlefont=[times,bold,18],
                color=lineColours,
                tickmarks= [ default,[seq( [j*(1+separation)=0, 1+j*(1+separation)=1][], j=0..columns-1)]],
                axis[1]= [gridlines= [[seq(j=j,j=0..rows)],colour= gray,linestyle= dash]]:
  timeValues:= seq( toPW(TT[j])+ (1+separation)*(j-1), j=1..columns,1):
  timePlots:= plot( ListTools:-Reverse([timeValues]), plotOptions):
  labelValues:= seq( [-1,0.5+j*(1+separation), yAxisLabels[j+1]], j=0..columns-1,1):
  labelPlots:= plots:-textplot(ListTools:-Reverse([labelValues]), font=[times,bold,12]):
  plots:-display(timePlots, labelPlots);

  separation:=1.5:
  plotOptions:= xAxis=0..rows-1, title = "Timing diagram",titlefont=[times,bold,18],
                color=lineColours,
                tickmarks= [ default,[seq( [j*(1+separation)=0, 1+j*(1+separation)=1][], j=0..columns-1)]],
                axis[1]= [gridlines= [[seq(j=j,j=0..rows)],colour= gray,linestyle= dash]]:
  timeValues:= seq( toPW(TT[j])+ (1+separation)*(j-1), j=1..columns,1):
  timePlots:= plot( ListTools:-Reverse([timeValues]), plotOptions):
  labelValues:= seq( [-1,0.5+j*(1+separation), yAxisLabels[j+1]], j=0..columns-1,1):
  labelPlots:= plots:-textplot(ListTools:-Reverse([labelValues]), font=[times,bold,12]):
  plots:-display(timePlots, labelPlots);
#
#  tabulate_id:= Tabulate(TT,width=40):
#  DocumentTools:-SetProperty([tabulate_id, 'fillcolor'[2 .. (), 2 .. ()], white],[tabulate_id,width,"40%"]):
#  DocumentTools:-SetProperty([tabulate_id, title,"Truth Table"],[tabulate_id,exterior,all]):
#  DocumentTools:-SetProperty([tabulate_id, fillcolor[1, 6], red],[tabulate_id,width,"40%"]):
#  DocumentTools:-SetProperty([tabulate_id, 'visible'[1 ..(), 1], false]):

 

 

 

 

 

 

 

Download logDiag4.mw

@Carl Love 

I did say there were many other ways to do it. Perhaps you find the attached more "acceptable"?

 

  restart;
  with(Logic):
  P:=TruthTable( X &xor Y &xor Z, form=MOD2):
  toPW:=M->piecewise(seq( [t<j, M[j]][], j=1..numelems(M))):
  plot( [ seq( toPW(P[j])+2*(j-1),
               j=1..4
             )
        ],
        t=0..8,
        tickmarks=[ default,
                    [seq( j=j mod 2, j=0..7)]
                  ]
      );

 

 


 

Download logDiag3.mw

First 85 86 87 88 89 90 91 Last Page 87 of 207