@emendes 

Please clarify the problem's size for the large example that you just ran. You referred to the "length" of the two input lists (which I'll call A and B throughout this thread). If this refers to Maple's length command, then that information is not very useful. Knowing the problem's "size" (via multiple metrics) is essential for guessing which algorithm will be most efficient.  So please execute the following, which should run quickly even for very large input.

ProbSize:= (A::list(set), B::list(set))->
    [nops~, (add@nops~)~, (nops@`union`@op)~]([A,B])
:
#For example, using your initially presented problem:
ProbSize(ans6, ans7); 
             [[20, 30], [120, 210], [16, 14]]

And if you still have access to the output of your test run (which I'll call L), please do the following:

S:= (add@nops~)(L); 
            S:= 47

The number S provides a theoretical lower bound on the time complexity of the computation.

Naive algorithms for this problem (such as yours, vv's, Tom's first, and mine above) have a runtime proportional to nops(A)*nops(B). (The storage required by vv's (but not the others) is also proportional to this, unfortunately.) For the example that you initially presented, this product is 600; for the example you just tested, it's about 6 billion (yikes!) if by "length" you meant nops.

To not be naive, let's consider the sparsity of the problem. The numbers that I just asked you to compute are attempts to measure the sparsity. Comparing 47 to 600, the example problem is quite sparse. Tom's second algorithm attempts to exploit the sparsity (whether he knew it or not). But his Maple implementation of the algorithm has a serious flaw that makes that procedure nearly worthless for large input. The flaw is only in the implementation, not the algorithm itself. That flaw is the loop 

for i in L1[j] do T[i]:=`union`(T[i],{j}) od

As is well known and often discussed here, it's extremely inefficient to construct a long list, sequence, or set by iteratively appending, the time to do so being proportional to the ​​​​​​square of the number of iterations.

The problem in Tom's code can be fixed, and here is my version that does that:

SubsetPairsSparse:= proc(A::list(set), B::list(set))
local 
    a, b, x, j, 
    T:= table([for x in `union`(A[],B[]) do x= table() od])
;
    for j,b in B do for x in b do T[x][j]:= () od od;
    T:= {indices}~(T, 'nolist');
    [for a in A do `intersect`(seq(T[x], x= a)) od]
end proc
:

Notice that neither this nor Tom's version of it do any looping process (for, seq, map, ~) whose number of steps is nops(A)*nops(B). Please try running that on a largish example, compare with my SubsetPairs1a, and report back. If the example that you use is not the one mentioned in your most-recent Reply, then please run my problem-sizie metrics on it also. If that goes well, we'll work on a threaded version.

Regarding the bipartite graph: Sorry, I forgot to relabel the vertices, which is needed for it to be recognized as bipartite. I've corrected the code above, and the entries from A will appear as negative numbers on their vertex labels. I also added to the DrawGraph command some basic options to improve readibility. Much more could be done along these lines.

Given an instantiation to numeric values as shown above, you may want to automatically verify whether those values are all real. That can be done like this (starting where the Answer above left off):

Inst:= eval(BoundV, FreeV) union FreeV;
type(Inst, set(name=realcons));
                              true

Can you provide a web page or paper that gives details of the algorithm? The information that I can find so far is too general to be of much use.

Just to clarify and to make sure: You're interested in an algorithm to factor multivariate polynomials with only rational coefficients, right?

Here is a complete Maple implementation covering all documented calling sequences of the Matlab function, and allowing Y to have any number (including 0) of dimensions, singleton dimensions, and starting indices other than 1. Of course, the vast majority of the code below is to handle all the weird cases; the mathematical computation of the discrete trapezoid rule is just a dot product.

trapz:= proc(x::rtable, y::{rtable, nonnegint}, dim::nonnegint, $)
description "Maple implementation of Matlab's trapz";
option `Author: Carl Love <carl.j.love@gmail.com> 2020-Dec-9`;
uses AT= ArrayTools;
local X, Y, Dim, nd, s, dims;
    #Deal with the 1-, 2-, or 3-argument calling sequences, and decide between
    #the 2 possible cases for the 2-argument sequence.
    if nargs=1 then Y:= x
    elif nargs=2 then 
        if y::rtable then (X,Y):= (x,y) else (Y,Dim):= (x,y) fi
    else (X,Y,Dim):= (x,y,dim)
    fi;

    #Get Y's dimensions:
    if (nd:= nops((dims:= [rtable_dims](Y))))=0 then return fi;
    local sz:= (rhs-lhs+1)~(dims);
    local LB:= (1-lhs)~(dims);

    #Find the key dimension of Y if it wasn't specified in the
    #arguments:
    if Dim::name then 
        if Y::Matrix then (Dim,s):= (2, sz[2])
        #Find 1st non-singleton dimension:
        else for Dim,s in sz while s=1 do od
        fi
    else
        try s:= sz[Dim]
        catch: error "Y has no non-singleton dimension"
        end try
    fi;

    #If X was passed in, check its size, and adjust its shape if
    #necessary; if not, set it to be unit spacing:
    if X::name then X:= `<|>`(1 $ (s-1))
    elif rtable_num_dims(X)<>1 or numelems(X)<>s then
        error "expecting X to be 1-D with %1 elements", s
    else 
        X:= AT:-Alias(X, [s], ':-row');
        X:= X[2..] - X[..-2]
    fi;

    #Prepare Y to be handled essentially as an Array of column vectors
    #(without actually breaking it up):
    local Yd:= AT:-Alias(Y/~2, sz, if nd=1 then ':-column' else fi);
    Yd:= 
        Yd[(()..()) $ (Dim-1), 2.., (()..()) $ (nd - Dim)] +~
        Yd[(()..()) $ (Dim-1), ..-2, (()..()) $ (nd - Dim)];
    s--;
    LB:= LB[[1..Dim-1, Dim+1..-1]];

    #final computation:
    if nd=1 then
        X.Yd
    else
        rtable(
            dims[[1..Dim-1, Dim+1..-1]][],
            (()-> 
                local J:= [args] +~ LB;
                X . AT:-Alias(
                    Yd[J[..Dim-1][], .., J[Dim..][]], [s], ':-column'
                )
            ),
            ':-subtype'= 
                if Y::Matrix then Vector[['-column', ':-row'][Dim]] 
                else Array
                fi
        )
    fi
end proc
:

@acer Okay, I understand your point. I was thrown off by the OP's use of the word precedence. Now that I reread the Question in light of your Reply, I see that the issue has nothing to do with precedence, and that the OP had simply misused that word (albeit in good faith).

@wswain You wrote:

• The original array elements were of a text(string) type, ...

Yes, that's correct.

• ... but you use subsindets to "transform" the list elements to string. 

No, that's incorrect. The second argument to subsindets, string in this case, tells the command to look for objects of type string in its first argument, the list. Hence, the list must already contain strings for this to work.

• Can you explain the "transformer" and "rest" arguments ...

The "rest" argument is not used in this example. If it had been used, it would be a constant (often an option) that would be the fourth argument of subsindets and passed as the second argument to the transformer.

• ... how you know a valid transformer operator (in this case 'string') ,,,

You mean operand, not "operator". The subsindets command only passes to the transformer objects of the correct type, string in this case, specified in its second argument.

• ...  vs it just acting on the elements?

Indeed, subsindets is overkill in this flat, linear case. It could be replaced (in Maple 2018 or later I believe) by

map(index, [seq](Array3[1, 2..]), [5..])

•  I presume a different substring operation could have been used to remove chars 1-4 ?

No, u-> u[5..] is precisely the substring operation that does remove characters 1-4.

The Answers by Tom and acer show that some relatively minor bug related to your specific Question was fixed by Maple 2016. However, the more-important bug, the operator precedence itself, was not fixed until Maple 2020 (I think at my request). Here's proof:

kernelopts(version);
  Maple 2019.2, X86 64 WINDOWS, Nov 26 2019, Build ID 1435526
a %* b %+ c %* d;
                     %*(%+(%*(a, b), c), d)

kernelopts(version);
  Maple 2020.1, X86 64 WINDOWS, Jul 30 2020, Build ID 1482634
a %* b %+ c %* d;
                     %+(%*(a, b), %*(c, d))

@Madhukesh J K At the moment, I think that you'll learn more by trying to follow written directions rather than an explicit code example. If you try this, and it doesn't work, please post an executed worksheet showing your attempt, and I'll give you a code example.

First, do you know what a coefficient is? And do you know what a nonlinear term in a differential equation is? If you're not sure about either of these, please ask. Otherwise, here are the explicit and simple steps:

1. Make C a coefficient of every nonlinear term in all 4 ODEs. In other words, multiply each such term by C.
2. Add option continuation= C to the dsolve command. 

That's all there is to it.

@Wilks You have taken a very easy problem and made it enormously complex. Like I said from the start, there's no need to use solve. Like I said after you provided details, there's no need to use int or intat. All that can be replaced by dsolve.

You need to put restart at the beginning of your worksheet.

The entire "chapter" of your worksheet entitled INTEGRACIÓ EQ. DE L'ELÀSTICA (LLEI DE DESPLAÇMENTS) can be replaced by these few lines of code, which computes all 8 Ys at once:

for Y in {Y||(1..8)} do
    unassign(Y); #just in case
    assign(
        Y= unapply(
            (rhs@dsolve)
                ({diff(Y(x), x$2) = cat(Y,`"`)(x), Y(0)=0, D(Y)(0)=0}),
            x
        )
    )
od:

That's if you want the Ys without the symbolic constants of integration. If you do want the constants, it can still be done with dsolve; let me know. It seems like you want the same Y's to have the constants in one part of the worksheet but not the other. Just use different names to distinguish them!

@David Sycamore Given any list L of [p,q] pairs (such as is generated by Kitonum's or my procedure), the sequence of just p can be extracted as Pseq:= op~(1, L)[]. The efficiency that can be gained by not storing the q in the first place is minimal.

@nm Time has been measured in minutes, so your last number should be "a little over 1/2 a minute". Otherwise, Vote Up.

@manju It wasn't me who posted that. However, that would be part of a command that sets all variables other than x and y in expression U to 1. The first part of the command would be eval(U, ...).

@mmcdara Yes, thanks for spotting that. I meant f(5). My error was the result of copying-and-pasting the old code and editing it (as is often the case with my posted errors). I'll correct the Reply above.

Type appliable is more general than type procedure, but it's often too-much-more general.

@Maple_lover1 Okay, my next guess is that the numbers that become strings are stored by Excel in some non-numeric format. If that is the case, the underlying issue cannot be addressed from the Maple side; however, it can be ameliorated by parse, as already shown.

@mmcdara Yes, I do prefer your "double-arrow procedure" method, and indeed I personally use that much more often than I use an indexed procedure.

So to apply this to the OP's situation:

div5:= proc(f) local x; unapply(f(x)/5, x) end proc:
f:= i-> x-> x^2+i:
div5(f(5));

Or, better yet, add some argument type checking, which isn't available for indexed procedures:

div5:= proc(f::appliable) local x; unapply(f(x)/5, x) end proc:
f:= (i::algebraic)-> (x::algebraic)-> x^2+i:
div5(f(5));

I see that the OP is using 2D Input. If instead they'd use 1D input (aka Maple Input), the 1st line could become (in Maple 2019 or later):

div5:= (f::appliable)-> local x; unapply(f(x)/5, x):