My own wish list starts with the large integer capabilities...
- A direct modular multiplicative inverse. Note that you can get exactly that from one of the returns of gca, however, a modular inverse makes a lot of sense, and is trivial to implement. (I've attached a simple modular inverse tool to this answer to show that it is indeed easy to compute, based purely on MATLAB's GCA utility.)
- A modular square root. This is not quite so trivial, but there are several fast algorithms to be found to compute, perhaps starting with Shanks-Tonelli.
- Modular n'th roots might be nice too, though it has been the sqrt that seems most valuable in my experience. I could survive without a modular n'th root.
- The fibonacci function is nice, but for some computations all you really need is mod(F_n,p), where p is some large integer. That is, while a student might like to see a million digit numberreturned, much of the time, you don't really need that million digit number returned, but just the modulus, and it is MUCH faster to compute the mod.
- Much like Fibonacci, the Lucas numbers can be valuable tools. They also should have the capability to be returned in a modular form.
- A direct quadratic residue computation would be nice, without needing to venture into the realm of MuPad.
Numbers 1 and 2 above can most easily be implemented as options to powermod. Thus we could allow a negative power for powermod, or a power of 1/2. We would have the modular inverse of x wrt base p as
powermod(x,-1,p)
Likewise, the modular square root of x would simply arise from
powermod(x,1/2,p)
The beauty of such an extension to powermod is it would require no new functions, just additional functionality for a code where it currently produces an error. powermod should and could give us everything we want.
Yes, I will concede that you can get a modular square root from mupad. For example...
X = 12345;
p = nextprime(1e9)
p =
1000000007
Xroot = feval(symengine, 'numlib::sqrtmodp', X, p)
Xroot =
258519153
mod(Xroot.^2,p)
ans =
12345
But I'm afraid this seems a bit of a kludge, requiring the use of feval. Anyway, the above computation will fail some of the time, seemingly returning a result when no modular root exists. And, yes, I could have predicted exactly that failure in advance.
X = 12346;
Xroot = feval(symengine, 'numlib::sqrtmodp', X, p)
Xroot =
419311443
mod(Xroot.^2,p)
ans =
999987661
feval(symengine, 'numlib::isquadres', X, p)
ans =
FALSE
with yet another call to the kludgistic feval.