Please see my blog post Gradients, Jacobians, Hessians, and higher order derivatives.

Inspired by Israel's post, define the inner product ip by:

ip := proc(f1,f2)
	option remember;
	if   type(f1,`+`) then map(ip,f1,f2)
	elif type(f2,`+`) then map2(ip,f1,f2)
	elif type(f1,`*`) then select(type,f1,freeof(t)) * ip(select(type,f1,dependent(t)),f2)
	elif type(f2,`*`) then select(type,f2,freeof(t)) * ip(f1,select(type,f2,dependent(t)))
	else int(int(f1 * f2,t = 0..tf),tf = 0..T)
	end if
end proc:

Then,

ip(f(t),g(t) + h(t));
ip(f(t) + g(t),h(t));
ip(2*sin(x)*f(t)*cos(y),g(t));
ip(f(t),2*sin(x)*g(t)*cos(y));

int(int(f(t)*g(t), t = 0 .. tf), tf = 0 .. T)+int(int(f(t)*h(t), t = 0 .. tf), tf = 0 .. T) int(int(f(t)*h(t), t = 0 .. tf), tf = 0 .. T)+int(int(g(t)*h(t), t = 0 .. tf), tf = 0 .. T) 2*sin(x)*cos(y)*int(int(f(t)*g(t), t = 0 .. tf), tf = 0 .. T) 2*sin(x)*cos(y)*int(int(f(t)*g(t), t = 0 .. tf), tf = 0 .. T)

Unless I am mistaken, your integral is proportional to the definite version of

int((a + b*cos(c/lambda)) / (1 + d*sin(e/lambda)^2),lambda);

where a, b, c, d, and e are constants (depending on k[2]). On face value, Maple is unable to integrate it, giving just int((a+b*cos(c/lambda))/(1+d*sin(e/lambda)^2), lambda) But maybe some hardcore mathematicians (I am only a physicist) can help you out here.

I do not know the exact reason for the error you get, besides that it seems to be related to the use of two levels of indexation, that is, k[n1][n2]. If reformulated in a matrix like way as

s:=k[1,1]*x[1]+k[2,2]*x[1];
s1:=k[1,1]*x[1]+k[1,2]*x[2];
l:=solve({coeff(s,x[1])=coeff(s1,x[1])},{k[1,1],k[2,2],k[1,2],k[2,1]});
p:=solve({coeff(s,x[2])=coeff(s1,x[2]),k[1,2]=0},{k[1,2],k[2,2],k[1,1],k[2,1]});
assign(p);

the problem seems to disappear. Is that useful for you?

I think the best you can do is to use the Physics package of Maple 11 as follows: First define a prefix for noncommutative quantities:

with(Physics):
Setup(noncommutativeprefix = N);

Then use the functions Commutator and Simplify as in the following examples:

• Containing a constant:

  Commutator(a * N1,N2);
  Commutator(N1,a * N2);
  Simplify(Commutator(a * N1,N2));
  Simplify(Commutator(N1,a * N2));
  

• Containing a sum:

  Commutator(N1 + N2,N3);
  Commutator(N1,N2 + N3);
  Simplify(Commutator(N1 + N2,N3));
  Simplify(Commutator(N1,N2 + N3));
  

• Containing a product:

  Commutator(N1 * N2,N3);
  Commutator(N1,N2 * N3);
  Simplify(Commutator(N1 * N2,N3));
  Simplify(Commutator(N1,N2 * N3));
  

Note that the quantities N1, N2, and N3 behave noncommutatively because of the prefix N, whereas a behaves commutatively (like a constant). Letting N1 and N2 denote your x and y, respectively, then your example may be correctly evaluated as

Commutator(c*N1 + d*N2,N2);

where c and d behave commutatively, as a above does.

Looking at the plot it is clear that the nontrivial intersection is (x,3) where x must be determined. As the graph of g goes linearly from (1,0) to (5,4), x must be given by (5-1)/(4-0)*3 + 1 = 4, using the formula x = 1/a*(y - y0) + x0, where a is the slope. PS: I did not know that piecewise functions could be entered the way you do - so thanks for that tip.

Using your example, let

A := {[1,0,1,0],[1,1,0,0],[0,0,1,1]}:

For any A, your requested equivalence classes are given by B below:

B := []:
for i from 1 to nops(A) do
for j from i+1 to nops(A) do
	s := A[i] + A[j];
	new := true;
	for k from 1 to nops(B) do
		if new and s = B[k][1] then
			B[k] := B[k] union {s};
			new := false
		end if
	end do;
	if new then B := [B[],{s}] end if
end do
end do:
B;

Or, much shorter, using the ListTools package:

with(ListTools):
L := [seq(seq(A[i] + A[j],j=i+1..nops(A)),i=1..nops(A))]:
B := map(convert,[Categorize((x,y) -> evalb({x[]} = {y[]}),L)],set);

The sign of lambda1 in the equation for N2 seems wrong: if the nuclei of type N1 decay to nuclei of type N2, then that sign should be positive.

ode := {
   diff(N1(t),t) = -lambda1*N1(t),
   diff(N2(t),t) = +lambda1*N1(t)-lambda2*N2(t)
}:
ics := {N1(0) = N,N2(0) = 0}:
sol := simplify(dsolve(ode union ics));

For some specific values for lambda1, lambda2, and N, the following plot may be obtained:

plot(subs({N = 1000,lambda1 = 2,lambda2 = 1},map(rhs,sol)),t = 0..5);

which seems physically correct.

I assume that the expression your have given should be the Faraday tensor of electromagnetism. However, in the present form, it is not; the signs of the two B1-terms must differ. I suppose, though, that it is a typo on your behalf, as is a missing );. With a corrected B1-sign, your Faraday tensor can be written as

with(DifferentialGeometry):
DGsetup([t,x,y,z],M):
F := evalDG(
   + E1*dt &w dx + E2*dt &w dy + E3*dt &w dz
   - B1*dy &w dz - B2*dz &w dx - B3*dx &w dy
);

where I have used the more powerful wedge-notation of differential calculus. The dual of any two-form might be defined as the map

Dual := (x) -> evalDG(
   + dx &w dy * Hook([D_z,D_t],x)
   + dy &w dz * Hook([D_x,D_t],x)
   + dz &w dx * Hook([D_y,D_t],x)
   + dt &w dx * Hook([D_y,D_z],x)
   + dt &w dy * Hook([D_z,D_x],x)
   + dt &w dz * Hook([D_x,D_y],x)
):

in which the two-form x is contracted, using the interior product Hook, with the dual basis (vector fields) D_t, D_x, D_y, and D_z. Using it on the Faraday tensor given, it yields your desired dual:

FDual := Dual(F);

I do not understand your basis. For the tensor you provide the basis must consist of three vectors. You must also provide either an explicit metric, or some inner product (for the space that the basis span) from which it can be obtained. The metric is needed because it is the inverted metric which must be used to raise the two indices of the tensor you provide.

Note added: Oh boy, I got carried away, only realizing after posting that my suggestion is nothing else than acer's function IR2 which I did not properly notice in the first place as my attention got entirely focused on his function IR. Here is another way:

M := (n::posint) -> Matrix(n,n,(i,j) -> `if`(i+j = n+1,-1,0)):
seq(M(i),i=1..5);

Concerning your first question, the appropiate function seems to be difforder:

with(PDEtools):
difforder(eq1);
                               1

Concerning your second question, the following function seems to solve the problem discussed in the thread Using select in conjunction with diff:

with(PDEtools):
findDiffOrder := (d,f) -> max(0,op(map(difforder,select(has,indets(d),f)))):
findDiffOrder(eq2,y),
findDiffOrder(eq2,u);
                              2, 1

Taking as an example the Taylor expansion of sin(x) about 0 to order 6 you can write something like

diff(taylor(sin(x),x=0,6),x$2);

Thanks for making it clear to me why my expectation is wrong. But then, please forgive me for saying so, I have a hard time figuring out how to use the packages in a truly unified way; it seems to me that I have to do my calculus at two places: in the Physics package if I want to do abstract index manipulation (by itself a very nice new feature, even though the ability to clearly distinguish between space and time is really missed), and in the old tensor package, for instance, if I want to do some 'number crushing'. Why this dichotomy? P.S.: This post should have been a reply to ecterrab, but became misplaced.

If you by decimal places mean the number of decimal places to the right of the decimal point, then the following function might be useful:

areEqual := proc(
	n::posint,
	p1::'set'([float,float,float]),
	p2::'set'([float,float,float])
)
	evalb(
		map(l -> map(f -> round(10^n*f),l),p1) =
		map(l -> map(f -> round(10^n*f),l),p2)
	)
end proc:

Note that it uses round. If you prefer you may use floor or ceil instead. The argument n is the number of decimal places to the right of the decimal point, and p1 and p2 are the two sets of lists of three floats. An example:

set1 := {[0.347,1.456,3.456],[2.2,3.3214,-4.521443]}:
set2 := {[0.347,1.456,3.456],[2.2,3.3216,-4.524673]}:
areEqual(1,set1,set2),
areEqual(2,set1,set2),
areEqual(3,set1,set2);
                       true, true, false