Carl Love

Carl Love

26832 Reputation

25 Badges

11 years, 289 days
Himself
Wayland, Massachusetts, United States
My name was formerly Carl Devore.

MaplePrimes Activity


These are answers submitted by Carl Love

The (kl)-shuffles that you mentioned in a Reply can be done by the oddly named command Iterator:-TopologicalSorts. The odd name comes from Donald Knuth; I don't know if it's standard. The command (when used with its inverse option) lists all permutations that satisfy certain positional inequalities. For example, the positional inequalities for a (2,3)-shuffle are {1 < 2, 3 < 4, 4 < 5}.

KLshuffle:= (K::posint, L::posint)->
local i, p;
    {seq}(
        [seq](p),
        p= Iterator:-TopologicalSorts(
            K+L, {seq(i < i+1, i= 1..K-1), seq(i < i+1, i= K+1..K+L-1)},
            'inverse'
        )
    )
:
KLshuffle(2,3);
     {[1, 2, 3, 4, 5], [1, 3, 2, 4, 5], [1, 4, 2, 3, 5], [1, 5, 2, 3, 4], 
      [2, 3, 1, 4, 5], [2, 4, 1, 3, 5], [2, 5, 1, 3, 4], [3, 4, 1, 2, 5],
      [3, 5, 1, 2, 4], [4, 5, 1, 2, 3]}

 

@MapleMathMatt Your procedure IsSubsequence is much faster, by a factor of about 3, than ArrayTools:-IsSubsequence. It can be simplified not making A into an Array and by not indexing it. Also, it's very easy to add an option to return the indices of such that B[K] = A when A as been verified as a sublist. (This is akin to the output= indices option of ArrayTools:-IsSubsequence.)   I did these things for code simplification and added functionality; the efficiency is about the same.

Definition: A list A is a sublist of a list B if there exists a strictly increasing list of positive integers (indices of B) such that A = B[K].

Here's my procedure:

Sublist?:= proc(A::list, B1::list, {indices::truefalse:= false}, $)
uses AT= ArrayTools;
local 
    n:= nops(B1), B:= rtable(B1), a, p:= 0,
    J:= [for a in A do if member(a, AT:-Alias(B, p, [p+1...n]), 'p') then p else return false fi od]
;
    if indices then J else true fi
end proc
:


Here's a time-test comparison of verifying a 10000-member sublist of a 10-million-member list:

R:= rand(1..10^4):
B:= ['R'()$10^7]:
K:= [combinat:-randcomb(10^7, 10^4)[]]:
A:= B[K]:

J:= CodeTools:-Usage(Sublist?(A, B, indices)):
memory used=79.58MiB, alloc change=76.30MiB, 
cpu time=78.00ms, real time=115.00ms, gc time=0ns

#Verify correctness:
type(J[..-2] - J[2..], list(negint)), evalb(B[J] = A);
                           true, true

CodeTools:-Usage(ArrayTools:-IsSubsequence(A, B, output= indices)):
memory used=114.70MiB, alloc change=-38.15MiB, 
cpu time=219.00ms, real time=471.00ms, gc time=93.75ms

CodeTools:-Usage(Sublist?(A, B)):
memory used=79.58MiB, alloc change=76.30MiB, 
cpu time=78.00ms, real time=104.00ms, gc time=0ns

CodeTools:-Usage(ArrayTools:-IsSubsequence(A, B)): 
memory used=114.58MiB, alloc change=114.45MiB, 
cpu time=312.00ms, real time=421.00ms, gc time=0ns

 

I strongly advise you to use a constructive approach (akin to what @sand15 showed) rather than a modification approach (akin to what @Kitonum showed). I am not considering any issues of computational effciency in this; the constructive approach leads to easier-to-read and easier-to-modify code.

If you were to use a modification approach, it can easily be done with subs without referring to any position within the list:

parms:= subs((thickness=0)= (thickness=3), parms);

However, I'd also want to avoid the specifc reference to 0. With that in mind, you can do

parms:= subsindets(parms, identical(thickness)=anything, ()-> thickness=5);

When options are repeated -- -such as your example that produces a list containing both thickness=0 and thickness=3 -- the option appearing at a higher position in the list takes precedence. In other words, there's no need to remove thickness=0 as long as thickness=3 appears anywhere in the list after it. This can be effectively utilized in plotting procedures by using the _rest meta-parameter. This lets the procedure's user specify ad hoc option modifications on the command line. For example:

MyPlot:= proc(MainArgs::list)
local default_opts:= (color= red, thickness= 0, linestyle= dash);
    plot(MainArgs[], default_opts,
_rest)
end proc
:
MyPlot([x^2, x= -2..2]);
MyPlot([x^2, x= -2..2],
thickness= 3);
 

Pad the text strings with the newline control character, which is "\n" in Maple. That will add a blank line either above or below the text snippet. And I agree with you totally that the relative linear measure of text should be determined by the font and be independent of the Cartesian coordinates of the base point of the text. Here's an example:

plots:-display(
    plot(
        [[[-3,1],[-1,1]], [[1,1],[3,1]], [[-3,-1],[-1,-1]], [[1,-1],[3,-1]]],
        thickness= 5
    ),
    plots:-textplot([
        [-2, 1, "A1", align= below],  [2, 1, "\nB2", align= below],
        [-2, -1, "C3", align= above], [2, -1, "D4\n", align= above]
    ]),
    axes= normal, size= [300$2]
);

There are many ways to represent a complex function in real 3D space. Here, I give separate animations of surface plots of the real and imaginary parts of your sol. Your t becomes the time-based animation parameter.

If you execute this code in a Maple worksheet, the plots will appear side-by-side with the t-values synchronized. It's difficult to copy that to this forum, so I present them here in separate Replies and not necessarily synchronized.

(plots:-display @ `<|>` @ plots:-animate~)(
    plot3d, 
    `[]`~(
        (evalc@[Re,Im])(exp(I*t)*sin(10*x)*sin(4*y)), x= 0..Pi, y= 0..Pi, 
        title=~ cat~(["Real", "Imaginary"], " part of solution\n")
    ),
    t= 0..1, frames= 49
);

real is a pre-defined Maple property; it doesn't become a property because you use it with assume or assuming. For example,

restart:
type(real, property);
              true

The `::` operator is an inert pair-forming operator. Just like the other such operatots (`=``<``..`, etc.), it takes exactly two operands of any type whatsoever (including NULL and multi-operand expression sequences), and evaluation of the expression is determined entirely by the surrounding context. Thus, the line you quoted from ?:: is wrong; indeed, that line is contradicted by the 7th paragraph in the same section of the very same help page:

  • In any other context, :: evaluates its arguments, but does not itself evaluate any further. Thus, it can be used as a data structure or part of a larger data structure. Although it is not required, it is customary to use :: as a data structure only when the right-hand side represents a type or a type-like concept. For example, see RealRange.

So, in `if`(a::b, c, d), the is expected to be a type. But in is(a::b), the b is first considered as a property. Many basic types can also be used as properties.

If you want to define a type real, it is easy to do so. I don't see any need for it, given the existence of type realcons. What did you have in mind?

Regarding the dictionary-style entries in Maple's help vs. Maple's actual help pages: Those entries give pure mathematical defintions, without regard to Maple or implementation; they're Platonic ideals (as it should be). Often, for practical reasons, the Maple implementation of a mathematical concept is not exactly the same as its pure definition. The same is likely true in any computer language. Floating-point number is a computer construct, not a mathematical one, so the pure-math dictionary has no entry for it.

@nm Types are syntactic predicates. They are not mathematical predicates. Unlike types, properties can be given to pure symbolic variables. For example, it makes no sense to say that x has type integer if x is an unassigned variable, because a symbolic variable is not an integer. On the other hand, x can be assume'd to have to the property integer and still remain a symbolic variable (integer is both a type and a property).

Here is the most important distinction between types and properties, and it's why there's no possibility of combining the concepts as you suggest. Types must be capable of returning either true or false after a finite number of computational steps. Returning FAIL is not a viable option. On the other hand, for a property of even trivial mathematical significance (such as being equal to 0), there will necessarily exist expressions with symbolic variables such that there is no computational means of determining whether the expression has the property. For these, is must return FAIL.

In summary:

  • Types are predicates to be applied to actual discrete structures existing in Maple's memory. Types should never return FAIL.
  • Properties are predicates about abstract, idealized (primarily pure math) structures. Properties must have the option of returning FAIL

As you said (using different words), TravelingSalesman returns a Hamiltonian cycle (returning to the starting point), whereas you want a Hamiltonian path (not returning). (Hamiltonian simply means that it includes every point.) A convenient way to do this is to add a dummy vertex (graph theory lingo for point) whose distance to all other vertices is a positive constant (say 1), then use TravelingSalesman, then discard the dummy. The fact that positioning such an equidistant dummy point is not physically possible in our ordinary conception of space is not significant. 

In the code below, n+1 is the dummy vertex. By using this as the startvertex and removing the first and last vertices from the returned cycle, we're left with a path containing just the original vertices. Note that the complete graph is constructed simply from its distance matrix (called a weight matrix in the lingo); nothing more elaborate than that is needed.

restart:
pts:= [[0,0], [1,2], [3,1], [5,4]]: 
n:= nops(pts):
GT:= GraphTheory: LA:= LinearAlgebra:
G:= (GT:-Graph@Matrix)(
   n+1, 
   (i,j)-> `if`(i=j, 0, `if`(i>n or j>n, 1, LA:-Norm(<pts[i]-pts[j]>, 2))), 
   shape= symmetric, datatype= hfloat  #This line is not essential
):
#[2..-2] removes 1st and last vertices from returned cycle:
HamPath:= GT:-TravelingSalesman(G, startvertex= n+1)[2][2..-2];

#The computation is done. All remaining code is just for plotting.
H:= GT:-InducedSubgraph(G, [$1..n]):  #Discard dummy
GT:-SetVertexPositions(H, pts);
GT:-HighlightTrail(H, HamPath);
GT:-DrawGraph(H, axes= frame, scaling= constrained);

There is an undocumented command for this: `tools/genglobal`.

restart:
sol:= x - x*`tools/genglobal`[1](_C)+x^2*`tools/genglobal`(_C);
new_constant:= `tools/genglobal`(_C);

 

Ronan: I think that the primary cause of your misunderstaning is not knowing what exactly is considered "invalid input". But before I explain that further, please read this explanation that I wrote for Preben regarding callseq_only:

@Preben Alsholm I think that you don't fully understand the purpose of callseq_only. (If you did already understand what I say below, then I apologize, and please let me know.) It's only purpose is to discriminate between different sources of errors beginning with "invalid input"; it has no influence on other errors. This is shown in the examples of help page ?overload.

The overload mechanism reacts only to errors that begin with "invalid input". The following 2 things are always true, regardless of any options:

  1. A procedure that is passed an argument that isn't allowed by its calling sequence will return an error beginning with "invalid input".
  2. A procedure may also return an error beginning with "invalid input" for some other reason. There are three possible ways (that I can think of right now): 
    a. The error came from another procedure called within the procedure.
    b. The programmer has explicitly put "invalid input ...." in an error command which was executed.
    c. (This is the case that is most relevant to this Question:) The error is caused by the procedure attempting to use a parameter for which no argument has been passed.

If a procedure has option overload, then we need to decide whether both 1 & 2 cause passage to the next procedure in the overload, or whether only 1 causes it. If you only want 1, then use option overload(callseq_only). Then type 2 errors (including type 2c) will be treated as normal errors, not triggers for the overload mechanism.

I don't think that there's any bug here, but I may change that opinion if you post your "contrived example".


Ronan: From reading the thread after Preben's Answer, I think that you now understand that failure to pass an argument to match a parameter is not necessarily an error (of any kind). And I think that you also now understand that it only becomes an error if the procedure's execution path tries to use the parameter whose argument is missing. Now, to complete your understanding, you need to know two more things:

  1. When such an error occurs, it is an "invalid input" error.
  2. Whether it triggers the overload mechanism depends on whether you used callseq_only.

Let me know if this does what you want:

`value/&||`:= (a,b)-> normal(a*b/(a+b))
:
`print/&||`:= proc(A,B)
uses T= Typesetting;
    T:-mrow(T:-Typeset(A), T:-mo(" &#8741; "), T:-Typeset(B))
end proc
:
R1 &|| R2 &|| R3;

     R1 || R2 || R3
value(%);

The &#8741; is the HTML code for the double vertical bar.

There is no list indexing needed at all:

restart;
X:= eval(
    <
        1,             -X3_3/2 - 1/2, 0,   -X2_3;
        -X3_3/2 - 1/2, -2*X3_4 - 1,   X2_3, 0   ;
        0,              X2_3,         X3_3, X3_4;
        -X2_3,          0,            X3_4,    1
    >,
    [X3_4, X3_3, X2_3]=~ [1-A^2, -1-A^2, A*(3*A^2 - 1)] /~ (3*A^2 - 1)
);     
andseq(
    andseq(
        is(e::nonnegative),
        e= LinearAlgebra:-Eigenvalues(eval(X, A= r))
    ), 
    r= solve(A^3 - A)
);

I wasn't sure about your question #3. The above returns true iff all 12 eigenvalues are real and nonnegative. Any variation of that that you want can also be done without indexing.

Note that is(e::nonnegative) is sufficient to check whether e is real and nonnegative.

Change dis play to display.

Here is a simple solution using select:

M3[select(M3[.., 3]=~ 2, [$op([1,1], M3)])];

To select columns whose 2nd row is 80 (for example), do

M3[.., select(M3[2]=~ 80, [$op([1,2], M3)])];

In either case, the = can be replaced by any of <<=>>=in, or <>.

For 2D Input, you need to put 1.. after $, as in 

M3[.., select(M3[2]=~ 80, [$1..op([1,2], M3)])];

Appending assuming complex is sufficient to remove the guaranteed nonzero factors in this case. A simple procedure for it is

Simp:= (q::`=`)-> local e:= factor((lhs-rhs)(q));
    `if`(e::`*`, (remove(f-> is(f<>0), e) assuming complex), e) = 0
:
eq:= diff(v(x),x,x)*exp(x^2) = 0:
Simp(eq);

You could change complex to real if you want.

Like this:

T:= Typesetting:
plots:-textplot(
    [1, 2, T:-mover(T:-msub(T:-mi("u"), T:-mn("1")), T:-mo("&rarr;"), 'mathcolor'= 'red')],
    'align'= {'above', 'right'}, 'font'= ["ARIAL", 'BOLD', 16]
);

3 4 5 6 7 8 9 Last Page 5 of 385