f=ones(31,1);
u=mgrid0(f);
