@J4James That would correspond to

plot3d(p(x)(y),x=0..1,y=-1.1..1.1,axes=boxed);

@Axel Vogt Your latest worksheet proves that the result is exactly zero, and that that fact was due to symmetries in a only.

I tried to simplify a little. The integrand considered as a function of 'a' really only depends on sin(2*a) and sin(a).

So the integrand can be written

u:=g(sin(2*a),sin(a));

where g is a function of two variables.
#Since sin is periodic with period 2*Pi we may integrate over -Pi..Pi instead of 0..2*Pi:
A:=Int(u,a=-Pi..Pi);
#Changing variables:
IntegrationTools:-Change(A,a=t-Pi);
#Changing the name of the variable:
subs(t=a,%);
#Changing interval again:
subsop(2=(a=-Pi..Pi),%);
#Thus 2*A is equal to
%+A;
#Which results in
   Int(g(sin(2*a), sin(a)), a = -Pi .. Pi)+Int(g(sin(2*a), -sin(a)), a = -Pi .. Pi)

#and that is zero if g is antisymmetric in its second variable, i.e. g(x,-y) = -g(x,y).

And our function in the concrete case satisfies that.

@Axel Vogt Your latest worksheet proves that the result is exactly zero, and that that fact was due to symmetries in a only.

I tried to simplify a little. The integrand considered as a function of 'a' really only depends on sin(2*a) and sin(a).

So the integrand can be written

u:=g(sin(2*a),sin(a));

where g is a function of two variables.
#Since sin is periodic with period 2*Pi we may integrate over -Pi..Pi instead of 0..2*Pi:
A:=Int(u,a=-Pi..Pi);
#Changing variables:
IntegrationTools:-Change(A,a=t-Pi);
#Changing the name of the variable:
subs(t=a,%);
#Changing interval again:
subsop(2=(a=-Pi..Pi),%);
#Thus 2*A is equal to
%+A;
#Which results in
   Int(g(sin(2*a), sin(a)), a = -Pi .. Pi)+Int(g(sin(2*a), -sin(a)), a = -Pi .. Pi)

#and that is zero if g is antisymmetric in its second variable, i.e. g(x,-y) = -g(x,y).

And our function in the concrete case satisfies that.

@J4James You are restricted to point style only if you plot M as is. And that would only be relevant if you imported M from a file. But then, as I said, you can recover the sequence dm, which was constructed directly from the plot data.

@J4James You are restricted to point style only if you plot M as is. And that would only be relevant if you imported M from a file. But then, as I said, you can recover the sequence dm, which was constructed directly from the plot data.

@Markiyan Hirnyk combine works too:

combine(t^n/t^k);

@Markiyan Hirnyk combine works too:

combine(t^n/t^k);

Could you elaborate?

Just one more:

evalindets(current_marks,list(integer),x->`+`(op(x))/nops(x));

Just one more:

evalindets(current_marks,list(integer),x->`+`(op(x))/nops(x));

@ilods It seems to be a simple zero/zero problem occurring when t (or s) is zero, so your "dirty" solution should be OK. However, the problem can be removed by inserting a line after the definition of k. Here is the line defining k followed by the new one:

k:=eval(B,{x=c[1],y=c[2]});
k:=map(normal~,k);

@ilods It seems to be a simple zero/zero problem occurring when t (or s) is zero, so your "dirty" solution should be OK. However, the problem can be removed by inserting a line after the definition of k. Here is the line defining k followed by the new one:

k:=eval(B,{x=c[1],y=c[2]});
k:=map(normal~,k);

@LijiH Yes, and it does so right away as you have also observed. But apparently you are not satisfied with that? Why not?

@LijiH Yes, and it does so right away as you have also observed. But apparently you are not satisfied with that? Why not?

@Chrissieh18 You could upload a worksheet.

Actually I strongly suspect that you are missing a multiplication sign between r and t in the exponential function. That is also the reason for the way I formulated my answer.