@epostma I'd probably suggest something like this, in case someone else is interested in the same:

module expand_units()
  local ModuleApply := proc(expr, $)
    return subsindets(expr, 'specfunc'({Unit, Units:-Unit}), worker);
  end proc;

  local worker := proc(u, $)
    local as_factors := convert(op(u), 'list', '`*`');
    local as_pairs := map(convert, as_factors, 'list', '`^`');
    return mul(map(pair -> Unit(pair[1]) ^ pair[2], as_pairs));
  end proc;
end module;

@C_R This is indeed a limitation in the system. You can see why below (the below is for Maple 2023 and beyond; the names were different in 2022 and before):

kernelopts(opaquemodules=false);
Units:-unitContextRecordTable["newton", "SI"]:-conversion;
# result: 1000*BaseUnitStruct(1)*BaseUnitStruct(2)/BaseUnitStruct(3)^2

The BaseUnitStruct index refers to the order in which base units were added to the system; Units uses Units:-backSubstitutionSet to turn it back into "regular" units. The representation simply does not allow associating the prefix with any particular base unit struct.

@dharr Unfortunately, this approach (and the shorter versions based on it, below) does not support unit annotations. (See e.g. https://www.maplesoft.com/support/help/maple/view.aspx?path=Units%2Fannotations.) If you use this procedure on Unit(mile / gal(petroleum)), it would attempt to turn petroleum and not gal into a separate unit, which doesn't work.

To fix this, you probably can't use subsindets in the inner transformer, which means that you'd need to actually write a procedure there. Not particularly difficult, especially using convert(..., list, `*`) and convert(..., list, `^`), but if I understand correctly that falls outside of C_R's original request.

@Nikol I'm sorry to hear you're still experiencing this problem. I wonder what's going on, though, because if I run the same code on my machine in Maple 2023, it works as expected. Here's a quick copy & paste from the TTY version:

$ maple2023 -B -s
    |\^/|     Maple 2023 (X86 64 LINUX)
._|\|   |/|_. Copyright (c) Maplesoft, a division of Waterloo Maple Inc. 2023 
 \  MAPLE  /  All rights reserved. Maple is a trademark of                 
 <____ ____>  Waterloo Maple Inc.                                          
      |       Type ? for help.         
> with(Statistics):
> Histogram(Sample(Poisson(200),10000));
memory used=4.7MB, alloc=41.3MB, time=0.19

0.03--------------------------***--------------------------------------+  
  +                           D DDD                                    |  
  +                         DDD D D                                    |  
  +                         D D D DDD                                  |  
0.025                       D D D D D                                  |  
  +                         D D D D D                                  |  
  +                      DDDD D D D D                                  |  
  +                      D  D D D D D                                  |  
0.02                     D  D D D D D  DDD                             |  
  +                      D  D D D D D  D D                             |  
  +                      D  D D D D DDDD D                             |  
  +                      D  D D D D D  D D                             |  
0.015                DDDDD  D D D D D  D DDD                           |  
  +                  D D D  D D D D D  D D D                           |  
  +                  D D D  D D D D D  D D D                           |  
0.01                 D D D  D D D D D  D D D                           |  
  +                DDD D D  D D D D D  D D DDD                         |  
  +                D D D D  D D D D D  D D D D                         |  
  +                D D D D  D D D D D  D D D DDDD                      |  
0.005           DDDD D D D  D D D D D  D D D D  D                      |  
  +           DDD  D D D D  D D D D D  D D D D  DDD                    |  
  +           D D  D D D D  D D D D D  D D D D  D DDD                  |  
  +       DDDDD D  D D D D  D D D D D  D D D D  D D DDDDD              |  
0 +********************************************************************+  
       160         180        200         220         240        260     
>                                                                        
 

The simplest way to do this is to apply an invertible linear transformation to your set of generators. In your example, you've taken x-z = (x-1) - (z-1), x+z-2 = (x-1) + (z-1), and y+z-2 = (y-1) + (z-1). That is, you've applied the linear transformation given by . Using any invertible linear transformation would give the same ideal. It's an interesting question how often this would transform a binomial set of generators into another binomial set of generators; if you're working over characteristic zero and with an ideal that cannot be generated by a single element, I'm pretty sure the answer is "almost always".

@C_R "display" is defined in the first section, which is folded in this MaplePrimes post but unfoldable in the Maple Cloud version here. (Or in the workbook version, downloadable from the top of this page or the Maple Cloud page.)

This is unfortunately not something we can support, because of the way that operator overloading works.

Here's what happens in detail. First of all, Maple's parser parses B - A as B + (-A). Then the Units:-Simple package's minus operator is evaluated on the result of evaluating A, that is, 1, 2, 3. Here we run into the first problem: we can't distinguish this call from a user calling `-`(1, 2, 3): we get the same argument sequence, there is no difference between the two. For the top level minus operator, this last call would be interpreted as 1 - 2 - 3, which evaluates to -4. So Units:-Simple's minus operator returns the same.

Now the plus operator gets called, with as its arguments the results of evaluating B and the earlier result that came out of the minus operator. B evaluates to 4, 5, 6, and the minus operator evaluated to -4. It is the nature of expression sequences in Maple that they automatically concatenate into a single, flat, expression sequence, so the kernel now passes this as the argument sequence 4, 5, 6, -4 to the plus operator. Again, we cannot distinguish this from a user calling `+`(4, 5, 6, -4). As the top level plus operator would return 4 + 5 + 6 + (-4) = 11 from that, so does the Units:-Simple package's plus operator.

Generally, arithmetic with expression sequences in Maple is extremely fragile. It is impossible for any overloaded package to replicate the top level operators exactly, and I don't believe there are documented rules for how the top level operators interpret expression sequences. Moreover, replacing the operator form with the functional form (that is, replacing something like A + B with something like `+`(A, B)) must also fail, because the expression sequence A, B immediately loses the information of where the split between the two was. This is the same issue we saw above: packages can really only provide the functional form of these operators.

I would strongly recommend against using expression sequences for arithmetic. Using lists or Vectors is much more reliable.

Erik Postma (he/him)
Manager, mathematical software group
Maplesoft

One thing I immediately noticed about the second listing is that, while correct, it is unnecessarily inefficient because it grows the list of primes one element at a time. That is quadratic in memory and time in Maple, because every new version of the list needs to be created from scratch. Instead it's better to use rtables (Arrays or Vectors) which can grow in-place. You can append with ArrayTools:-Append or the ",=" operator.

I wonder what ChatGPT would say to that...

@Carl Love @Angelo Melino  Thank you for finding this! I just submitted a fix for the main development version of Maple; it should make it into the next major version (which might be released next spring).

Erik Postma
Manager, mathematical software group

@Anthrazit The 'eval' should be done when reading from the data table; when you write it, the unit expressions will hopefully be valid -- otherwise you would get errors at that point, anyway, whenever you do anything with them -- so then the eval is not necessary.

@vv Thank you for this request for high performance, high precision evaluation of the Zeta function. That'll help us prioritize where we focus our efforts. If we hear this from more people, that would be a strong indicator to help redirect effort to this project.

Erik Postma
Manager, mathematical software group.

@Axel Vogt - one doesn't, in this case. It was me who suggested this to the support team, and I thought the `in` operator was defined in the PolynomialIdeals[Operator] subpackage, but it's actually in the PolynomialIdeals package itself.

Yes, that works for the Bootstrap command, thanks for pointing that out!

The focus of the user's question was on generating samples, so I answered that question initially. I believe they wanted to roll their own bootstrap-like mechanism. For just generating (large) samples, EmpiricalDistribution works well enough, and efficiently.

I figured I could show that you can use Maple's Bootstrap command, too, and because I had already constructed the EmpiricalDistribution object, that's where my mind went, and I didn't think about the fact that Bootstrap accepts the plain data sample, too.

By the way -- I think we can probably speed Bootstrap up a lot by having it sample into a larger buffer. It currently just uses the sample size for the buffer, so it fills a size-4 buffer 10000 times. It would be much more efficient to fill a size-10000 buffer 4 times, and then slice it into 2500 samples with ArrayTools:-Alias. I expect you would get performance very close to the calling sequence  you show in that case.

Thanks Robert, I enjoyed reading this account.

Two notes, one mathematical and one Maple-related. First the mathematical one: I think you would expect 3 degrees of freedom, right? There are 10 equations, but there is a dependency between them: all four row-sums must needs sum to the same as all four column-sums. So there are really only 9 independent equations.

The Maple point is that the combinatorial question you raise is solved by the Iterator package, in particular the BoundedComposition command. You supply it with a k-tuple of bounds (b1, ..., bk), and a sum s, all nonnegative integers. It then finds all k-tuples (a1, ..., ak) of nonnegative integers with ai <= bi that sum to s. This raises one problem: we want integers from 1-9, and  BoundedComposition will also give us zeroes -- but that's easily dealt with: we just use integers from 0-8 instead, adjust the requested sum (by subtracting k=4), and presto. For example:

with(Iterator):
k := 4;
B := BoundedComposition([8 $ k], 25-k):
Number(B); # result: 324
Print(B, 'showrank');

tells us that there are 324 such tuples with sum 25, and will print all of them. To do a computation with them, you would use:

for b in BoundedComposition([8 $ k], 25 - k) do
  tup := b +~ 1;
  # Now use the tuple, tup, which will contain suitable integers 1-9. (It's an Array.)
end do;

Of course, if one number is given already, you could generate 3-tuples.

@C_R If this is what you want, the best option is probably to use convert/system to convert every unit to one consistent system, and then use Carl's solution below:

expr := sin(t * Unit(1/min)) * Unit(mm);
to_SI := evalindets(expr, specfunc(Units:-Unit), e -> convert(e, 'system', 'SI'));
result := eval(to_SI, Units:-Unit = 1);
# result is now sin(t/60) / 1000