What were you expecting the entries of the completed matrix to look like? Obviously they can't be numbers since you've not defined W or S. But if they're not numbers, what's the point of your evalf? Are you trying to set up some finite-difference equations?

To other readers: Yeah, I know that the proximate cause of the error is the failure to evaluate the derivatives before substituting numbers for the variables. I don't think that this OP could make use of that information before the questions above are addressed.

@digerdiga You sound as if you're upset that something that most users would consider a bug has been fixed. Do you not like it? I think that almost all elementary calculus instructors would consider the old behavior a bug. Note that conversion in the other direction (which seems more important to me) still works the same way:

convert(1/(1-x), FormalPowerSeries);

@mmcdara 

Here's a sleeker and more-optimized version of yesterday's module, still done with just arithmetic and list operations. I put this into your time-testing worksheet and it beat all the other methods. IIRC, it's a factor of 7 or 8 faster than my module from yesterday. Please test/confirm that.

NimMatrix:= module()
local
   P, N, NK, 
   ModuleApply:= proc(N::And(posint, Non(1)), K::posint)
   local i,j;
      thismodule:-N:= N;
      P:= N^~[$0..K-1];
      NK:= P[-1]*N;
      `mod`:= modp;
      baseN:= proc(n) option remember; local q:= n; ['irem(q, N, 'q')'$K] end proc;              
      rtable(
         'symmetric', (1..NK)$2, 
         [seq([seq(i &Nim j, i= 0..j)], j= 0..NK-1)], 
         'subtype'= 'Matrix'
      )
   end proc
;
export baseN, `&Nim`:= (a,b)-> inner(baseN(a)+baseN(b) mod N, P);             
end module:

The mysterious command inner in procedure `&Nim` is an undocumented built-in command that computes the inner product (aka dot product) of two lists as if they were vectors.

Why specifically do you want 4 iterations?

Your first integral is utterly trivial. Your second integral is very close to this simpler one (they only differ by a constant factor):

I3:= int(sin(n*x)*csc(x), x);

I haven't been able to get Maple to do that one. If someone else can get Maple to do it, nearly the same steps should work on your integral. If n is any specific integer, then Maple can do it. So, we just need a summation formula for sin(n*x) assuming n::integer.

@mmcdara 

You wrote:

• Right, it was only to say that we are not doing all the operations in Z/pZ (p being the base) but that some steps are done in an another base (in fact the base in which the elements of the matrices are written)

The correct operations are done in Z/pZ digit-wise (or element-wise). There is no "another base" that is relevant to this problem. You are confusing an abstract positive integer---which is independent of human culture and has no "base"---with a system for writing those integers---which is culturally dependent.

@mmcdara There you go again: saying that my answer "wasn't correct" and that I "had corrected an error"! And, no, it is not due to vagueness on your part about the carry rule. Tom was using my original code, which is based only on the algorithm stated in the 4 steps in the second paragraph of your original Question. The problem is that there's an error in your statement of those 4 steps.

I must emphasize that my original code wasn't correct due to some kind of default behavior. On the contrary, it's clear from the way that I sectioned the code that I expended extra effort to write it to perform the steps as you wrote them, even though I had some skepticism that that was how you intended them. That skepticism was due to mathematical intuition telling me that there's no serious mathematical significance to base 10, but it was ameliorated by the fact that people often post puzzle-type problems that do rely on a base-10 representation.

And don't go thinking that I'm one of those people who never admits to having made an error. The record here on MaplePrimes shows numerous times that I have admitted to making an error without the slightest bit of argument about it.

@mmcdara I did not "correct" anything! My original procedures correctly performed the operations that you incorrectly stated in the Question. It was only after my "probing" that I was able to deduce what the operations were that you actually intended. This was only a deduction of intentions; you still haven't anywhere in this thread made a complete and correct statement of the algorithm implied by your Question, although you do have some code that correctly performs it.

Your Question mentions "base 10" three times. This problem has nothing whatsoever to do with base 10. You'll notice that nowhere in any of the correct code is the number 10 used. The sentence that you used that leads to the most confusion is "do c = a+b just as if a and b were numbers in base 10." Think about what that sentence means when the target base is greater than 4 until you understand that my original code correctly performs that operation. See also my second Reply to Tom below.

Thank you for the timing code; I was curious.

@mmcdara But you do have the foreknowledge that I sought: For any given matrix, K is fixed, and every list or vector used during the construction of the matrix can be length K. That's much more efficient than zero-padding each list separately.

The base being prime has nothing to do with the construction of these matrices using my method. Whether you choose to only construct them when the base is prime is up to you.

So, here's my optimized code:

NimMatrix:= module()
local
   sav_tablesize,
   ModuleLoad:= proc() sav_tablesize:= interface('rtablesize'); return end proc,
   ModuleUnload:= proc() interface('rtablesize'= sav_tablesize); return end proc,

   ModuleApply:= proc(N::And(posint, Non(1)), K::posint)
   uses AA= ArrayTools:-Alias;
   local 
      P:= AA(Array(0..K-1, k-> N^k), 0, [K]),
      NK:= P[-1]*N #NK = N^K
   ;
      interface('rtablesize'= max(interface('rtablesize'), NK));
      #The procedures baseN and `&Nim` are defined here so that 
      #they get fresh remember tables for each Matrix.
      baseN:= 
         proc(n::nonnegint)
         option remember;
         local q:= n, k;
            <seq(irem(q, N, 'q'), k= 1..K)>            
         end proc
      ;     
      `&Nim`:= 
         proc(a::nonnegint, b::nonnegint)
         option remember;
            local r:= add(irem~(add(baseN~([a,b])), N)*~P);
            #Implement symmetry:
            if a=b then r else procname(b,a):= r fi
         end proc
      ;         
      AA(Array((0..NK-1)$2, `&Nim`), 0, [NK$2])
   end proc
;
export
   baseN, `&Nim`
;
   ModuleLoad()
end module:
 

Usage:

M:= NimMatrix(3,3);

@tomleslie The reason for the different answers between our procedures is that my procedure Nim is doing the list addition "with carry". I did this because the OP stated incorrectly in the Question that the addition of the digits was to be done in base 10. It's actually not supposed to be done in any base, just natural integer arithmetic. The base-10 part seemed so weird to me that I got a clarification about "carry" from the OP. With the clarification, I now know that there's no carry, and that This algorithm has absolutely nothing to do with base 10! Some people confuse the natural unbased integer arithmetic with base 10. It only superficially seems that way because the end results are converted to base 10 for display purposes only. Otherwise, n is not the same thing as convert(n, base, 10), obviously.

This issue only matters, and our procedures will only differ, when the base is greater than 4.

I have a new batch of more-efficient procedures that I'm about to put up.

@mmcdara If you want to go the polynomial route, use PolynomialTools:-FromCoefficientList rather than series.

My B10 would be faster as a for loop because then I'd not need to compute each power of b separately. In my Answer, I was going for brevity of code (without sacrificing any generality).

If you have foreknowledge of the overall maximum list length needed, then ZPad can be made much faster by using Vectors of that length. It appears that for your work that there is such foreknowledge.

I need for you to clarify if a "carry" operation for the addition is needed when the base is greater than 4. Suppose that the base is 9 and the original operands are 7 and 8. Do we do this:
[[7], [8]] #operands converted to base 9
[15]       #sum of the above
[6]         #above mod 9
6           #above reconverted to base 10

or this:
[[7,0], [8,0]] #operands converted to base 9, padded for carry
[5,1]            #sum of the above
[5,1]            #above mod 9
14               #above reconverted to base 10
?

As soon as you can clarify that, I'll be able to give you an efficient procedure to create a whole Matrix.

Hmm. I just realized that I could figure it out directly from your code, presuming that you think that your code works for the higher-base cases. And it appears that you use the first of my two options. In that case, my Nim procedure won't work for bases greater than 4, and is generally inefficient.

@tomleslie The operation a +~ b will fail if the two lists have different lengths. See procedure ZPad in my Answer.

@tomleslie The OP is developing an algorithm for a BVP solution method that's not currently implemented in Maple. Of course, that algorithm is ultimately intended to be used for problems for which there's no easy symbolic solution. Obviously this work is in its early stages. A good way to test it is to compare its results against an exact result that's known to be correct.

This is all me guessing based on the OP's questions over the past few days.

@acer I mean that the execution of the OP's error statement prevents my procedure's return value from being displayed because now the OP's wrapping procedure has no return value.

Is that clear now? Sorry, I often convey the context by using boldface, which may not be apparent in a literal or out-loud reading. When I use boldface here on MaplePrimes, it's always meant to represent literal Maple code. If that Maple code in also an English word that's meaningful in the same context, I can see that that can cause confusion. So, I should've said "blocked by the execution of an error statement" rather than "blocked by error."

@VEIND You've quickly come to a situation requiring unapply. It's due to the order in which things are evaluated. In this code:

f:= x-> diff(g(x), x);
f(0);

the order is not as it appears on screen. The order is

1. Make the unevaluated expression diff(g(x), x) into a function of x, and assign it to f
2. Replace x with 0
3. Evaluate the resulting expression diff(g(0), 0), which is nonsense because you can't differentiate with respect to a number.

Changing that to

f:= unapply(diff(g(x), x), x);
f(0);

makes the order

1. Evaluate diff(g(x), x), and make the result a function of x assigned to f
2. Replace x with 0, and do the arithmetic.