@Markiyan Hirnyk 

First, I build list  L . It is a list of frequencies for the sum of the digits in a sequence of six digits,  ie how many times this or that sum of of digits occurs.

L = [[0, 1], [1, 6], [2, 21], [3, 56], [4, 126], [5, 252], [6, 462], [7, 792], [8, 1287], [9, 2002], [10, 2997], [11, 4332], [12, 6062], [13, 8232], [14, 10872], [15, 13992], [16, 17577], [17, 21582], [18, 25927], [19, 30492], [20, 35127], [21, 39662], [22, 43917], [23, 47712], [24, 50877], [25, 53262], [26, 54747], [27, 55252], [28, 54747], [29, 53262], [30, 50877], [31, 47712], [32, 43917], [33, 39662], [34, 35127], [35, 30492], [36, 25927], [37, 21582], [38, 17577], [39, 13992], [40, 10872], [41, 8232], [42, 6062], [43, 4332], [44, 2997], [45, 2002], [46, 1287], [47, 792], [48, 462], [49, 252], [50, 126], [51, 56], [52, 21], [53, 6], [54, 1]]

For example, the sub-list  [30, 50877]  means that the sum of digits equal to 30 in the 6-digit sequence occurs 50877 times. The same means the list  M , only for 5 digits.

Note that 23 = 6 +6 +6 +5. Therefore  I find all representations of  N , which is given as argument of the procedure, as a sum of 4 numbers (the first number from 0 to 54, and so on). Get the list  S .

Then I go through the list  S in the cycle , and, for example for N = 100, element of the list  [20, 35, 30, 15]  in accordance with the well-known rule of combinatorics (rule of product) gives me  35127*30492*50877*13992  variants.

 .

@Markiyan Hirnyk 

First, I build list  L . It is a list of frequencies for the sum of the digits in a sequence of six digits,  ie how many times this or that sum of of digits occurs.

L = [[0, 1], [1, 6], [2, 21], [3, 56], [4, 126], [5, 252], [6, 462], [7, 792], [8, 1287], [9, 2002], [10, 2997], [11, 4332], [12, 6062], [13, 8232], [14, 10872], [15, 13992], [16, 17577], [17, 21582], [18, 25927], [19, 30492], [20, 35127], [21, 39662], [22, 43917], [23, 47712], [24, 50877], [25, 53262], [26, 54747], [27, 55252], [28, 54747], [29, 53262], [30, 50877], [31, 47712], [32, 43917], [33, 39662], [34, 35127], [35, 30492], [36, 25927], [37, 21582], [38, 17577], [39, 13992], [40, 10872], [41, 8232], [42, 6062], [43, 4332], [44, 2997], [45, 2002], [46, 1287], [47, 792], [48, 462], [49, 252], [50, 126], [51, 56], [52, 21], [53, 6], [54, 1]]

For example, the sub-list  [30, 50877]  means that the sum of digits equal to 30 in the 6-digit sequence occurs 50877 times. The same means the list  M , only for 5 digits.

Note that 23 = 6 +6 +6 +5. Therefore  I find all representations of  N , which is given as argument of the procedure, as a sum of 4 numbers (the first number from 0 to 54, and so on). Get the list  S .

Then I go through the list  S in the cycle , and, for example for N = 100, element of the list  [20, 35, 30, 15]  in accordance with the well-known rule of combinatorics (rule of product) gives me  35127*30492*50877*13992  variants.

 .

@Markiyan Hirnyk 

Your example is not valid because:

1) Each numeric multiset should be consisting of 6 numbers.

2) Both multisets should be equal.

No one taught me to work with Maple. This applies also to writing programs and procedures in Maple. 

This is not my work, but just hobby in spare time. My profession - is math and I teach the students 

on this subject.
I think that you can not learn to write programs just by reading books. 

Of course, you need to learn some simple rules, but most important you should to solve specific problems,

and refer to the books as needed.

Just discovered this topic and I noticed that the code in my original post again distorted. I would love to know why this happens? The correct code is copied into the comment so that everyone can compare with the code from the original post. Also mw-file was uploaded.

Polygons:=proc(N, S)

local Area, n, m, It, rk, L, U, T0, f, T1;

global T;

uses ListTools;

if type(N, odd) or S<(N-2)/2 or S>N^2/16 then

print(`No solutions`) else

Area:=proc(L::list)

local M;

M:=[op(L), L[1]];

add(M[i,1]*M[i+1,2]-M[i+1,1]*M[i,2], i=1..nops(L))/2;

end proc:

n:=(N-2)/2; m:=floor((n+1)/2);

It:=proc(K)

local j, k, r, M, s;

j:=0; k:=nops(K[1]);

for r in K do

if r[k]=[0,0] then M:=[[0,1], [1,0]] elif

r[k]=[0,m] then M:=[[0,m-1], [1,m]] elif

r[k]=[n,0] then M:=[[n-1,0], [n,1]] elif

r[k] in {seq([i, n+1-i], i=m..n)} then M:=[[r[k,1]-1,r[k,2]], [r[k,1],r[k,2]-1]]

elif

r[k,1]=0 and r[k,2]>0 and r[k,2]<m then M:=[[0,r[k,2]-1], [0,r[k,2]+1], [1,r[k,2]]] elif

r[k,2]=m and r[k,1]>0 and r[k,1]<n+1-m then M:=[[r[k,1]-1,m], [r[k,1],m-1], [r[k,1]+1,m]] elif

r[k,2]=0 and r[k,1]>0 and r[k,1]<n then M:=[[r[k,1]-1,0], [r[k,1],1], [r[k,1]+1,0]] else

M:=[[r[k,1]-1,r[k,2]], [r[k,1],r[k,2]+1], [r[k,1]+1,r[k,2]], [r[k,1],r[k,2]-1]]; fi;

rk := r[1..k-1];

j := j+1;

s[j] := seq(`if`(member(p,rk), NULL, [op(r),p]), p = M);

od;

[seq(s[i], i=1..j)];

end proc;

L:=[seq([[0,i]], i=0..m)];

U:=(It@@(N-1))(L);

T0:=[seq(`if`(((p[1,1]=p[N,1] and abs(p[N,2]-p[1,2])=1) or (p[1,2]=p[N,2] and abs(p[N,1]-p[1,1])=1)) and abs(Area(p))=S, p, NULL),  p = U)];

f:=(x,y)->convert(x, set)=convert(y, set);

T1:=[Categorize(f, T0)];

T:=[seq(T1[i,1], i=1..nops(T1))]: # List of all the polygons

nops(T);  # Number of all the polygons

fi;

end proc:

Polygons_of_matche.mw

@Carl Love

I don't know the reason of this. In the morning when I edited my original post, it was the same thing, but in the evening, when I opened this page everything looked normal. 

@Carl Love

Thanks for the discovery of an analogy with the structures of polyominoes. Of course, if we have a description of all polyominoes of given type (fixed area) in terms of the coordinates of the vertices of the corresponding polygon, we can calculate it's areas through the coordinates of the vertices (by Area procedure) and to select the desired polygons.
I usually place the pictures on the website radical.ru and then use green tree in the corner to give the link.

@Markiyan Hirnyk

You are right. The main idea of ​​the solution is the same as in the topic  http://www.mapleprimes.com/questions/142492-How-Many-Routes-In-A-4x4-Square#comment142513  . Very useful formula for calculating the area of ​​a simple polygon (not self-intersecting) by the coordinates of its vertices  ( Area  procedure). In addition to reducing the amount of computations, I ultimately narrowed down the search area. It is easy to prove that any simple polygon bounded  N  matches can be placed in a polygon with vertices  (0, 0),  (n, 0),  (n, 1),  (m, n),  (0, m) , where  n=(N-2)/2  and  m=floor((n+1)/2) .

@Markiyan Hirnyk 

Thanks for the marked discrepancy. This happens when there is only one variant as in your example  N=10, S=6 .

Replace the lines

T1:=Categorize(f, T0);

T:=[seq(T1[i,1], i=1..nops([T1]))]:

by

T1:=[Categorize(f, T0)];

T:=[seq(T1[i,1], i=1..nops(T1))]:

 Now the procedure is working correctly in all cases.  As for visualization, the code for it is an particular case, and it  is written for N=12, S=6 only (here there are 35 variants).

@Markiyan Hirnyk 

In the example, replace  Polygon  into  Polygons . I found no link to edit my post. 

PS. Link found and the bug was fixed. 

for the detected error and improving of the algorithm! After correcting my code, I got the same result. Now there are no duplicates:

Routes(4, 6);
nops({op(T)});

      68272

      68272

for the detected error and improving of the algorithm! After correcting my code, I got the same result. Now there are no duplicates:

Routes(4, 6);
nops({op(T)});

      68272

      68272

Check the specification of your equation!

You also did not answer these questions: "What means your dual equality? The system of equations?"

Check the specification of your equation!

You also did not answer these questions: "What means your dual equality? The system of equations?"

Your code does not work, if the list  r  contains several  elements with the same maximum frequency:

r := [3, 3, 3, 3, 4, 4, 4, 3, 3, 3, 3, 3, 3, 3, 3, 2, 4, 3, 1$13];

with(Statistics):

T := Tally(r);

lhs(op(select(c-> is(rhs(c) = max(map(rhs, T))), T)));

       r := [3, 3, 3, 3, 4, 4, 4, 3, 3, 3, 3, 3, 3, 3, 3, 2, 4, 3, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1]

                                            T := [1 = 13, 2 = 1, 3 = 13, 4 = 4]

Error, invalid input: lhs expects 1 argument, but received 2

Your code does not work, if the list  r  contains several  elements with the same maximum frequency:

r := [3, 3, 3, 3, 4, 4, 4, 3, 3, 3, 3, 3, 3, 3, 3, 2, 4, 3, 1$13];

with(Statistics):

T := Tally(r);

lhs(op(select(c-> is(rhs(c) = max(map(rhs, T))), T)));

       r := [3, 3, 3, 3, 4, 4, 4, 3, 3, 3, 3, 3, 3, 3, 3, 2, 4, 3, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1]

                                            T := [1 = 13, 2 = 1, 3 = 13, 4 = 4]

Error, invalid input: lhs expects 1 argument, but received 2