get_matrix_fem.m

五点差分型多重网格方法：各种插值算子的比较）

%== nx,ny is numbers of elements in two directions respectively =========
%
function [x,y,Ag,N,bg]=get_matrix_fem(nx,ny)
N = (nx+1)*(ny+1);
gauss = [-1/sqrt(3), 1/sqrt(3)]; % Gaussian quadrature point
% construct Q1 mesh
[x,y,conn,ne,np] = getQ1mesh(nx,ny);
Ag = zeros(np);
bg = zeros(np,1);
nloc = 4; % number of nodes per element


for ie = 1:ne % loop over all elements
    
   [A_loc,rhse] = elemA(conn,x,y,gauss,ie);
   % assemble local matrices into the global matrix
   for i=1:nloc
       irow = conn(ie,i); % global row index
       bg(irow)=bg(irow) + rhse(i);
       for j=1:nloc
           icol = conn(ie,j); %global column index
           Ag(irow, icol) = Ag(irow, icol) + A_loc(i,j);
       end
   end
end
Ag = sparse(Ag);