matlab在[0,1]x[0,1]上生成六角形区域

根据点，把六条线求出来，然后生成点对，最后根据线性约束把想要的点提取出来 

clear all
clc
close all
N = 3000;
x = rand(1, N);
y = rand(1, N);
xy = [x;y];
k1 = 1;
domain = zeros(2, N);
for i =1:Nif (xy(2,i)>0.25) && (1.5*xy(1,i)+0.25-xy(2,i)>0) && (-1.5*xy(1,i)+1.75-xy(2,i)>0)domain(1,k1)=xy(1,i);domain(2,k1)=xy(2,i);k1=k1+1;end
end
k2=k1;
for i =1:Nif (xy(2,i)<0.75) && (-1.5*xy(1,i)+0.75-xy(2,i)<0) && (1.5*xy(1,i)-0.75-xy(2,i)<0)domain(1,k2)=xy(1,i);domain(2,k2)=xy(2,i);k2=k2+1;end
endinteriork=1;
for kk=1:Nif domain(1,kk)~=0 && domain(2,kk)~=0interiorD(1,interiork)=domain(1,kk);interiorD(2,interiork)=domain(2,kk);interiork=interiork+1;end
end%增加边界点
M = 50;
xbt1 = linspace(0,1.0/3, M);
ybt1 = ones(1,M)*0.75;
kbt=k2;
for i=1:Mdomain(1,kbt)=xbt1(i);domain(2,kbt)=ybt1(i);kbt=kbt+1;
endxbt2 = linspace(1.0/3,0.5, M);
ybt2 = 1.5*xbt2 + 0.25;
for i=1:Mdomain(1,kbt)=xbt2(i);domain(2,kbt)=ybt2(i);kbt=kbt+1;
endxbt3 = linspace(0.5,2.0/3, M);
ybt3 = (-1.5)*xbt3 + 1.75;
for i=1:Mdomain(1,kbt)=xbt3(i);domain(2,kbt)=ybt3(i);kbt=kbt+1;
endxbr1 = linspace(2.0/3,1.0, M);
ybr1 = ones(1,M)*0.75;
kbr=kbt;
for i=1:Mdomain(1,kbr)=xbr1(i);domain(2,kbr)=ybr1(i);kbr=kbr+1;
endxbr2 = linspace(5.0/6,1.0, M);
ybr2 = 1.5*xbr2 -0.75;
for i=1:Mdomain(1,kbr)=xbr2(i);domain(2,kbr)=ybr2(i);kbr=kbr+1;
endxbr3 = linspace(5.0/6,1.0, M);
ybr3 = (-1.5)*xbr3 + 1.75;
for i=1:Mdomain(1,kbr)=xbr3(i);domain(2,kbr)=ybr3(i);kbr=kbr+1;
endxbb1 = linspace(2.0/3,1.0, M);
ybb1 = ones(1,M)*0.25;
kbb=kbr;
for i=1:Mdomain(1,kbb)=xbb1(i);domain(2,kbb)=ybb1(i);kbb=kbb+1;
endxbb2 = linspace(0.5,2.0/3, M);
ybb2 = 1.5*xbb2 -0.75;
for i=1:Mdomain(1,kbb)=xbb2(i);domain(2,kbb)=ybb2(i);kbb=kbb+1;
endxbb3 = linspace(1.0/3,0.5, M);
ybb3 = (-1.5)*xbb3 + 0.75;
for i=1:Mdomain(1,kbb)=xbb3(i);domain(2,kbb)=ybb3(i);kbb=kbb+1;
endxbb1 = linspace(2.0/3,1.0, M);
ybb1 = ones(1,M)*0.25;
kbb=kbr;
for i=1:Mdomain(1,kbb)=xbb1(i);domain(2,kbb)=ybb1(i);kbb=kbb+1;
endxbb2 = linspace(0.5,2.0/3, M);
ybb2 = 1.5*xbb2 -0.75;
for i=1:Mdomain(1,kbb)=xbb2(i);domain(2,kbb)=ybb2(i);kbb=kbb+1;
endxbb3 = linspace(1.0/3,0.5, M);
ybb3 = (-1.5)*xbb3 + 0.75;
for i=1:Mdomain(1,kbb)=xbb3(i);domain(2,kbb)=ybb3(i);kbb=kbb+1;
endxbL1 = linspace(0,1.0/3, M);
ybL1 = ones(1,M)*0.25;
kbL=kbb;
for i=1:Mdomain(1,kbL)=xbL1(i);domain(2,kbL)=ybL1(i);kbL=kbL+1;
endxbL2 = linspace(0.0,1.0/6, M);
ybL2 = 1.5*xbL2 + 0.25;
for i=1:Mdomain(1,kbL)=xbL2(i);domain(2,kbL)=ybL2(i);kbL=kbL+1;
endxbL3 = linspace(0.0,1.0/6, M);
ybL3 = (-1.5)*xbL3 + 0.75;
for i=1:Mdomain(1,kbL)=xbL3(i);domain(2,kbL)=ybL3(i);kbL=kbL+1;
endirk=1;
for kk=1:Nif domain(1,kk)~=0 || domain(2,kk)~=0irregularD(1,irk)=domain(1,kk);irregularD(2,irk)=domain(2,kk);irk=irk+1;end
endfigure('name','domian')
scatter(irregularD(1,:),irregularD(2,:),'r.')
hold onfigure('name','u')
U = sin(pi*irregularD(1,:)).*sin(pi*irregularD(2,:));
size2U = size(U);
c = U;
scatter3(irregularD(1,:),irregularD(2,:),U,20,c,'.');
grid on
colorbarXY = irregularD;