If you have an account on old mapleprimes.com, login there and click on Maple 9.5 Add-on Software Components.

listplot3d is not the command you need - it should be pointplot3d. For example,

plots[pointplot3d](readdata("DATA",3), axes=boxed, labels=[x, y, z],
connect=true);

For example,

f:=(k1,k2,n1,n2)->add( add( `if` (j1*k1/n1 +j2*k2/n2+k2/n2 > 1/2 and
j1*k1/n1 +j2*k2/n2 >1/2, 1, 0), j1=0..n1), j2=0..n2):

f(1,2,3,4);

                                  17

evalf uses an approximate value of Pi that produces a relatively small imaginary component in sin. To get rid of it, one can apply fnormal to the result and use simplify/zero after that.

For example, if n=b then print("n may not be completely factored") end if; or if n=b then printf("%d may not be completely factored ",n) end if; depending on what you want to see printed. Also, userinfo can be used.

What do you expect to see on the plot? It has some sense to plot integrals of them though, i.e., sums of Heaviside functions. For example,

a:=inttrans[fourier](sin(x),x,s);

               a := Pi (Dirac(s + 1) - Dirac(s - 1)) I

plot(int(a/I,s),s=-2..2);

You provided the data for xi and yi, but forgot to include the values f(xi,yi) that should be fitted by P(x,y). Well, here is an example of fitting values of i,

X:=<150, 152, 154, 156, 156, 157, 158, 158>:
Y:=<134, 136, 106, 61, 90, 124, 41, 135>:
Z:=<($0..7)>:
Statistics:-Fit(a*x+b*y+c*x*y+d,<X|Y>,Z,[x,y]);

  1.9657693669022 x + 1.30452326136629115 y

         - 0.00821594396071781072 x y - 304.913905308599738

The simplest thing, probably, is adding few print statements in your procedure - for example, time() can be printed in front of loops that you suspect may cause problems, and after them (or maybe few times inside) together with some labels (loop names). Also, userinfo can be used, profile or debugger. See corresponding help pages, ?userinfo, ?profile and ?debugger. See also ?CodeTools,Profiling.

for i to 100 while not i=5 or for i to 100 while i<>5 Be aware though that both of these commands do the loop only 4 times - for i equal to 1, 2, 3, 4, and stop after the while condition is not satisfied, a:=0: for i to 100 while not i=5 do a:=a+i od: a; 10 If one wants to continue the loop after that, then instead of using while, if should be used inside a loop, such as in the following example, a:=0: for i to 100 do if i<>5 then a:=a+i fi od: a; 5045

Just skip style. Also, S can be defined more simple using seq, S:=[seq([i,irem((i+3)^2,7)],i=0..150,5)]; plot(S);

Let's try to find that experimentally,

m:=Statistics:-NonlinearFit(a*x,[1,2],[3,7],x,weights=[1,10],output=solutionmodule):
f:=unapply(m:-Results(leastsquaresfunction),x):
m:-Results(residualsumofsquares);

                           0.24937655860348

(f(1)-3)^2+100*(f(2)-7)^2;

                             0.2493765586

It looks as if the least-square error is defined as

                         Sum(w[i]^2*(f(x[i])-y[i])^2,i)

Test it on another example,

m:=Statistics:-NonlinearFit(a*x^2+b*x,[2,3,5],[3,7,12],x,weights=[2,5,-1],
output=solutionmodule):
f:=unapply(m:-Results(leastsquaresfunction),x):
m:-Results(residualsumofsquares);

                           7.7575757575758

evalhf(4*(f(2)-3)^2+25*(f(3)-7)^2+(f(5)-12)^2);

                         7.75757575757575158

Correct!

I would rather use something like that,

Abs:=rcurry(LinearAlgebra:-VectorNorm,2):

If abs is preferable, it can be done using a module,

Abs:=module() export abs; 
abs:=v->if v::Vector then LinearAlgebra:-VectorNorm(v,2) 
else :-abs(args) fi end:

Now, after with(Abs); the command abs will calculate the Euclidean norm for Vectors and work as usual otherwise. Similarly other Maple commands can be replaced - doesn't matter built-in or not.

That, probably, should go to the Everything but... forum. Anyway, here is a link.

You can start from reading ?liesymm and ?tensor. There are some other useful packages, such as GRTensor and Clifford.

?fnormal says that one should use simplify(expr,zero) for that.