function [A, B, C, num, bound] = NaD(gx, gy, NOE)  %NaD函數用來劃分網格gx與gy分別為橫縱方向的網格數,NOE為三角元總數
nx = gx + 1;
ny = gy + 1; %nx,ny分別為橫縱方向的結點數
num(ny, nx) = 0; %用來存儲節點號的數組
if(ny > nx)%用來存儲邊界節點的數組
    bound(4, ny) = 0; 
else
    bound(4, nx) = 0;
end
%將各個結點的編號存入num數組中
for j = 1 : nx
    for i = 1 : ny
        if (j == 1)
            num(i, j) = i; %第一列的結點編號
        else
            num(i, j) = num(i, j - 1) + ny; %非第一列的節點編號
        end
    end
end

%將邊界結點的編號存入bound數組中，其中1，2，3，4行分別對應為上邊界下邊界左邊界和右邊界
for j = 1 : nx
    bound(1, j) = num(1, j); %上邊界
    bound(2, j) = num(ny, j);%下邊界    
end

for i = 2 : ny - 1
    bound(3, i) = num(i, 1); %左邊界
    bound(4, i) = num(i, nx);%右邊界
end

A(1) = 1; B(1) = 2; C(1) = ny + 2; A(2) = 1; B(2) = ny + 2; C(2) = ny + 1;
epc = 2 * gy;
%計算奇數號三角元的頂點編號
for e = 3 :2 : NOE
    for e1 = 3 : 2 :e
        if (mod(A(e1 - 2), ny) == ny - 1)
            A(e1) = A(e1 - 2) + 2;
            B(e1) = B(e1 - 2) + 2;
            C(e1) = C(e1 - 2) + 2;
        else
            A(e1) = A(e1 - 2) + 1;
            B(e1) = B(e1 - 2) + 1;
            C(e1) = C(e1 - 2) + 1;
        end
    end
end
%計算偶數號三角元的頂點編號
for e = 4 :2 : NOE
    for e1 = 4 : 2 :e
        if (mod(A(e1 - 2), ny) == ny - 1)
            A(e1) = A(e1 - 2) + 2;
            B(e1) = B(e1 - 2) + 2;
            C(e1) = C(e1 - 2) + 2;
        else
            A(e1) = A(e1 - 2) + 1;
            B(e1) = B(e1 - 2) + 1;
            C(e1) = C(e1 - 2) + 1;
        end
    end
end
end