Use the  eliminate  command for this. The reverse representation is ambiguous and depends on the choice of parameters. For your example, if we take  x  and  y  as parameters, it is trivial:

Download elim.mw

A:=<1,2; 3,4>;
B:=<5,6; 7,8; 9,10; 11,12>;
<A[1], B[2], B[4]>;

An alternative to the  map  command would be the  seq  command:

restart;
L:=[[1,2,3],[7,8,9],[13,12,11]]:
[seq(l[3], l=L)];
L1:=[1, [2, 3], [4, [5, 6], 7], [8, 3], 9]:
[seq(`if`(l::list,l[1],l), l=L1)];

Addition. You do not need to form a system of equations in any form. Using the  LinearAlgebra:-LinearSolve  command you can solve the system simply by specifying the matrix of coefficients in front of the unknowns and the vector of the right side:

M:=(Matrix([[rho*h*omega[m,n]^2+k[1, 1],k[1,2],k[1,3],k[1,4]],[k[2,1],rho*h*omega[m,n]^2+k[2,2],k[2,3],k[2,4]],[k[3,1],k[3,2],k[3,3]+rho*h*omega[m,n]^2,k[3,4]],[k[4,1],k[4,2],k[4,3],k[4,4]+rho*h*omega[m,n]^2]]));
V:=-E*(Vector(1..4,[k[1,5],k[2, 5],k[3,5],k[4,5]]));
LinearAlgebra:-LinearSolve(M, V);

Download solve_(1)_new.mw

The  AlmostCarmichaelNumbers  procedure finds all the numbers that meet your requirements for a number  k  and for a range  R  (which you specify when using) .

AlmostCarmichaelNumbers:=proc(k::posint,R::range)
local T, N, m, n, L;
T:=table(); N:=0; m:=lhs(R);
for n from `if`(m::odd,m,m+1) by 2 to rhs(R) do
L:=ifactors(n);
if nops(L[2])=k then if `and`(map(t->is(t[2]=1),L[2])[]) and `+`(seq(`if`(irem(n-1,t[1]-1)=0,1,0),t=L[2]))=k-1 then N:=N+1; T[N]:=n fi; fi;
od;
convert(T, list);
end proc:

Examples of use for k=3, 4, 5 :

AlmostCarmichaelNumbers(3,1..10000);
AlmostCarmichaelNumbers(4,1..10000);
AlmostCarmichaelNumbers(5,1..10000);
              

A:=(-108*q - 12*sqrt(12*p^3 + 81*q^2))^(1/3)/6;
B:=(A^3)^(1/3);
c:=op([1,2], B);
subsop([1,2]=-(c^2)^(1/2), B);

Addition. This method can be used in some cases when simplifying radicals, for example:

restart;
A:=sqrt(11-4*sqrt(6))+sqrt(21+6*sqrt(6));
sqrt(A^2);

The differential equation must be entered with explicit indication of the derivative (not differentials) using Maple syntax. For your example, you can do this:

restart;
de:=2*y*dx+3*x*dy=0;
subs([dy=diff(y(x),x)*dx,y=y(x)], de);
dsolve(%, y(x));

Just differentiate your equation 2 times, taking  y as a function of  x :

restart;
Curves:=x^2+A*y(x)^2=B:
Eq:=diff(Curves,x,x);
dsolve({Eq,y(0)=sqrt(B/A),D(y)(0)=0}, y(x), implicit);
solve(%, B)=B;

Edit. Sure enough the right solution is Rouben's one, not me. I forgot to eliminate the constant  A . Rouben very politely "did not notice" this error.

You can do this (like as in the original post) without  DataFrame  , using only matrices:

A:=Matrix(3, (i,j)->i);
DF:=< <<D1|D2|D3>, A> | <``,D1,D2,D3> >;

You can use the  Explore  command for this:

Explore(plot([sqrt(a*x^3+b*x^2+c*x+d),-sqrt(a*x^3+b*x^2+c*x+d)], x=-3..3, color=red, view=-7..7), a=-3...3., b=-3...3., c=-3...3., d=-3...3.);

Maple's possibilities are not limitless, so just reduce the ranges for calculations. It is easy to check that one of the roots is exact 0:

Download calculs_new.mw

You can use the inert form  Diff  that will not be calculated, but can be calculated at any time with the command  value :

u[1](T[1], T[2], T[3]):=exp(T[1]*I):
Diff(diff(u[1](T[1], T[2], T[3]), T[1]), T[2]);
value(%);

Of course, you can make inert differentiation for both variables:

u[1](T[1], T[2], T[3]):=exp(T[1]*I):
Diff(u[1](T[1], T[2], T[3]), T[1], T[2]);
value(%);

Here are 2 customised legends with vertical alignment (text rotation is available only from Maple 2018):

restart;
with(plots):
Plots:=plot([t^2,t^3], t=-2..2, color=[blue,red], linestyle=[1,3], view=-2..2, labels=[t,"x,y"]):
Lines:=plot([[1.6,t,t=0.9..1.3],[1.9,t,t=0.9..1.3]], color=[blue,red], linestyle=[1,3]):
Labels:=textplot([[1.55,1.55,x=t^2,rotation=Pi/2],[1.85,1.55,y=t^3,rotation=Pi/2]]):
display(Plots, Lines, Labels);

Use  seq  command for this.

Example:

seq(sin(x), x=0..evalf(Pi/4), 0.1);

     0, 0.9983341665e-1, .1986693308, .2955202067, .3894183423, .4794255386, .5646424734, .6442176872

Check:=proc(n)
local L;
L:=ifactors(n);
if n>1 and type(n,odd) and (not isprime(n)) and `and`(map(t->is(t[2]=1),L[2])[]) and `or`(map(t->is(irem(n-1,t[1]-1)=0),L[2])[]) then true else false fi;
end proc: 

Example of use:
Check(561);
                             true

The previous answer skipped checking of square free.