f := int((A*sin(omega*t+phi)-B*sin(omega*t))^2, t):

g := int((A*sin(omega*t)*cos(phi)+A*cos(omega*t)*sin(phi)-B*sin(omega*t))^2, t):

delta=simplify(expand(f)-expand(g));  # Difference between f and g 

We see that the difference between  f  and  g  does not depend on  t, so there is no contradiction.

ex := taylor(exp(x), x = 0, 5);

eval(ex, x=1);

subs(x=1, convert(ex, polynom));

A:=Matrix(3,2, symbol=a):

B:=Matrix(2,3, symbol=b):

C:=Matrix([[8, 2, -2], [2, 5, 4], [-2, 4, 5]]):

simplify(B.A, {seq(seq((A.B)[i,j]=C[i,j], j=1..3), i=1..3)});

mahmood180, I tried to upload your file, but the 404 error appeared.

Look at my solution:

restart;

SpecialMatrix := proc (A::list)

local r, a, b, A1, A2, a1, a2; 

r := (1/2)*nops(A)-1/2;

assign(seq(a[i-1] = A[i], i = 1 .. r+1), seq(b[i] = A[r+1+i], i = 1 .. r));

A1[1, 1] := Matrix(1, {(1, 1) = 2*a[0]});

A1[1, 2] := Matrix([[seq(a[i], i = 1 .. r)]]);

A1[1, 3] := Matrix([[seq(b[i], i = 1 .. r)]]); A1[2, 1] := A1[1, 2]^%T;

A1[2, 2] := Matrix(r, {seq(seq((i, j) = a[j-i], j = i+1 .. r), i = 1 .. r-1), seq((i, i) = 2*a[0], i = 1 .. r)}, shape = symmetric);

A1[2, 3] := Matrix(r, {seq(seq((i, j) = b[j-i], j = i+1 .. r), i = 1 .. r-1), seq((i, i) = 0, i = 1 .. r)}, shape = antisymmetric);

A1[3, 1] := A1[1, 3]^%T; A1[3, 2] := -A1[2, 3]; A1[3, 3] := A1[2, 2];

a1 := `<,>`(seq(`<|>`(A1[i, 1], A1[i, 2], A1[i, 3]), i = 1 .. 3)); 

A2[1, 1] := Matrix(1);

A2[1, 2] := A1[1, 2];

A2[1, 3] := A1[1, 3];

A2[2, 1] := Matrix(r, 1);

A2[2, 2] := Matrix(r,{seq(seq((i, j) = a[i+j], j = 1 .. r-i), i = 1 .. r-1)});

A2[2, 3] := Matrix(r,{seq(seq((i, j) = b[i+j], j = 1 .. r-i), i = 1 .. r-1)});

A2[3, 1] := A2[2, 1]; A2[3, 2] := A2[2, 3]; A2[3, 3] := -A2[2, 2];

a2 := `<,>`(seq(`<|>`(A2[i, 1], A2[i, 2], A2[i, 3]), i = 1 .. 3)); 

(1/2)*a1+(1/2)*a2; 

end proc:

Example for r=4:

SpecialMatrix([seq(a[i], i = 0 .. 4), seq(a[i]^`&lowast;`, i = 1 .. 4)]);

SymmetricSum:=proc(S::{list,set}, P::list)

local L;

uses combinat;

if nops(S)<>nops(P) then error "Should be nops(S)=nops(P)" fi;

L:=permute(P);

add(mul(S[i]^l[i], i=1..nops(S)), l=L);

end proc;

CyclicSum:=proc(S::list, P::list)

local n, L0, L;

if nops(S)<>nops(P) then error "Should be nops(S)=nops(P)" fi;

n:=nops(P);

L0:=[op(P),op(P)]; L:=[seq([seq(L0[i+j], j=0..n-1)],i=1..n)];

add(mul(S[i]^l[i], i=1..nops(S)), l=L);

end proc:

Examples:

SymmetricSum([x,y,z], [2,3,5]);

CyclicSum([x,y,z,u], [2,3,4,5]);

Formal parameters:  L - the list of your functions  psi ,  m - a symbol or a number,  M - a positive integer.

MMatrix := proc (L::list, m, M::posint)

Matrix(M, [seq([seq(L[i]((1/2)*(2*j-1)/m), j = 1 .. M)], i = 1 .. M)]);

end proc:

Example:

MMatrix([seq(psi[1, i], i = 0 .. 4)], 10, 5);

Of course instead of arbitrary functions psi, you can write as the entries of the list L your specific functions.

Another way of filling the interior of the cardioid used. It is based on an approximation by polygons. The same method was used in my procedure Picture. See  http://www.mapleprimes.com/posts/145922-Perimeter-Area-And-Visualization-Of-A-Plane-Figure-

restart;

N := 120:

r := 2-2*sin(2*Pi*k/N):

L := [[0, 0], seq([r*cos(2*Pi*k/N), r*sin(2*Pi*k/N)], k = 0 .. N)]:

A := seq(plottools[polygon](L[1 .. n], color = cyan), n = 2 .. N+2):

B := seq(plottools[curve](L[2 .. n], color = black, thickness = 3), n = 2 .. N+2):

l := plottools[line]([0, 0], [2, 0], color = cyan, thickness = 4):

C := seq(plots[display](l, A[i], B[i]), i = 1 .. N+1):

plots[display](C, insequence = true, axiscoordinates = polar, scaling = constrained);

tmp:=[[0, 1, 2], [1, 0, 2], [1, 1, 2], [1, 2, 0]];

select(x->is(`+`(op(x))=3), tmp);

                            tmp := [[0, 1, 2], [1, 0, 2], [1, 1, 2], [1, 2, 0]]

                                       [[0, 1, 2], [1, 0, 2], [1, 2, 0]]

triType:=proc(a,b,c)

local A, B, C;

if is(a+b<=c) or is(a+c<=b) or is(c+b<=a) then error "No triangle" fi;

A:=b^2+c^2-a^2; B:=a^2+c^2-b^2; C:=a^2+b^2-c^2;

if is(A>0) and is(B>0) and is(C>0) then return acute else

if is(A=0) or is(B=0) or is(C=0) then return right else

obtuse fi; fi;

end proc;

Examples:

triType(3, 3, 4);

triType(3, 4, 5);

triType(3, 4, 6);

triType(3, 5, 10);

                                  acute

                                   right

                                 obtuse

Error, (in triType) No triangle

restart;

B:= ln(exp(t)+a+sqrt(exp(t)-a+b)):

e:=solve(A(t)=B, exp(t)):

subs(exp(t)=e[1], diff(B, t));

simplify(x^4+y^2+2*y*x^2, {x^2+y=z});

                                  z^2

or

algsubs(x^2+y=z, x^4+y^2+2*y*x^2);

                                  z^2

or

compoly(x^4+y^2+2*y*x^2, {x, y});   # composition of polynomials

                        x^2, x = x^2+y

For numerical solutions see  ?plots[odeplot]  and  ?DEtools[DEplot] . If you got the solution explicitly, ie. {x(t)=..., y(t)=..., ...} then  assign(%)  and then  plot([x(t), y(t)],...)  .

I do not think that  PlotPositionVector can solve your problem. A special procedure is required.
The procedure  CurvVector  returns  the list of coordinates the normalized curvature vector for a curve defined by the list of coordinates  L :

restart;

CurvVector:=proc(L)

local e1, e1n, e2, e2n, t;

uses Student[VectorCalculus], LinearAlgebra;

t:=indets(L)[1];

e1:=diff(L,t):

e1n:=unapply(simplify(expand(e1/Norm(convert(e1, Vector),2))), t) assuming t::realcons;

e2:=expand(diff(L,t,t)-(convert(diff(L,t,t),Vector).convert(e1n(t), Vector))*e1n(t));

e2n:=unapply(simplify(expand(e2/Norm(convert(e2, Vector),2))), t) assuming t::realcons;

end proc;

Your example:

with(Student[VectorCalculus]):

L:=[p*cos(p^2), p*sin(p^2)]:

f:=CurvVector(L);

V:= RootedVector(root=eval(L, p=1.5), -f(1.5));

A:=PlotVector(V, color=green):

B:=plot([p*cos(p^2), p*sin(p^2), p=1..2], thickness=2):

plots[display](A, B, scaling=constrained);

restart;

f:=x->tan(x)-x-1:  g:=D(f):

x[0]:=1.5: x[1]:=x[0]-f(x[0])/g(x[0]): epsilon:=10^(-3):

for n while evalf(abs(x[n]-x[n-1]))>=epsilon do

x[n+1]:=x[n]-f(x[n])/g(x[n])

od:

n;  # number of iterations

x[n];  # the root

                                                    7

                                          1.132267726

Verification:

fsolve(tan(x)-x-1);

                                           1.132267725

Actual accuracy is 10^(-9)

restart;

P := x->x^7+14*x^4+35*x^3+14*x^2+7*x+88;

x0 := (1+sqrt(2))^(1/7)-(-1+sqrt(2))^(1/7)-(3+2*sqrt(2))^(1/7)-(3-2*sqrt(2))^(1/7);

is(P(x0) = 0);

factor(P(x), (1+sqrt(2))^(1/7)):

x := solve(op(3, %));

is(x = x0);