function [A, b] = Assembler_1(p, e, t);
%clc
%clear all
%[p,e,t]=S_shapeadaptiv(2);

nnodes = size(p, 2);  % Antal noder. (= antal kolonner i p)
nedges = size(e,1);  % Antal triangelsidor som ligger på ett randsegment.
ntri = size(t, 2);    % Antal trianglar. (= antal kolonner i t)
%A = zeros(nnodes, nnodes);      % Initiera styvhetsmatris.
A=sparse(nnodes,nnodes);
b = zeros(nnodes, 1);           % Initiera lastvektor.
rv = zeros(nnodes, 1);          % Initiera randvektor.
q=0;
for n = 1:ntri

% Ta reda på hörnens nodnummer:  
  
  i = t(1, n);          % Nodnummer för hörn nummer 1 i triangel nummer n.
  j = t(2, n);          % Nodnummer för hörn nummer 2 i triangel nummer n.
  k = t(3, n);          % Nodnummer för hörn nummer 3 i triangel nummer n.
                        %
		        % 1 -> 2 -> 3 moturs riktning.
  
% Ta reda på hörnens koordinater:		       
		       
  xnodi = p(1, i);      % x-koordinat för nod nummer i.
  xnodj = p(1, j);      % x-koordinat för nod nummer j.
  xnodk = p(1, k);      % x-koordinat för nod nummer k.
  
  ynodi = p(2, i);      % y-koordinat för nod nummer i.
  ynodj = p(2, j);      % y-koordinat för nod nummer j.
  ynodk = p(2, k);      % y-koordinat för nod nummer k
 
 % Räkna ut triangelns area:
 
triangelarea = 0.5.*abs((xnodj - xnodi)*(ynodk - ynodi) - (xnodk - xnodi)*(ynodj - ynodi));
Triangelarea(n) = 0.5.*abs((xnodj - xnodi)*(ynodk - ynodi) - (xnodk - xnodi)*(ynodj - ynodi));
              
  
% Räkna ut (de konstanta) gradienterna av phi_i, phi_j, och phi_k:
  Dphi_iDx = (ynodj-ynodk)./-((xnodj-xnodk).*(ynodi-ynodk)-(xnodi-xnodk).*(ynodj-ynodk));          %(ynodj - ynodk)/(2*triangelarea);
  Dphi_iDy = (xnodk-xnodj)./-((xnodj-xnodk).*(ynodi-ynodk)-(xnodi-xnodk).*(ynodj-ynodk));          %(xnodk - xnodj)/(2*triangelarea);
  gradphi_i = [Dphi_iDx; Dphi_iDy];
  
  Dphi_jDx = (ynodk-ynodi)./-((xnodk-xnodi).*(ynodj-ynodi)-(xnodj-xnodi).*(ynodk-ynodi));          %%(ynodk - ynodi)/(2*triangelarea);
  Dphi_jDy = (xnodi-xnodk)./-((xnodk-xnodi).*(ynodj-ynodi)-(xnodj-xnodi).*(ynodk-ynodi));          %%(xnodi - xnodk)/(2*triangelarea);
  gradphi_j = [Dphi_jDx; Dphi_jDy];
  
  Dphi_kDx = (ynodi-ynodj)./-((xnodi-xnodj).*(ynodk-ynodj)-(xnodk-xnodj).*(ynodi-ynodj));   %(ynodi - ynodj)/(2*triangelarea);
  Dphi_kDy = (xnodj-xnodi)./-((xnodi-xnodj).*(ynodk-ynodj)-(xnodk-xnodj).*(ynodi-ynodj));    %(xnodj - xnodi)/(2*triangelarea);
  gradphi_k = [Dphi_kDx; Dphi_kDy];
  
  
  
  
  anodi =1; 
  anodj =1; 
  anodk =1; 

  fnodi = f(xnodi,ynodi);32.*xnodi.*(1-xnodi)+32.*ynodi.*(1-ynodi);% Värdet i nod i av f.
  fnodj = f(xnodj,ynodj);32.*xnodj.*(1-xnodj)+32.*ynodj.*(1-ynodj);  % Värdet i nod j av f.
  fnodk = f(xnodk,ynodk);32.*xnodk.*(1-xnodk)+32.*ynodk.*(1-ynodk);  % Värdet i nod k av f.
 
 A(i,i)=A(i,i)+gradphi_i'*gradphi_i*triangelarea;%./3;
 A(i,j)=A(i,j)+gradphi_i'*gradphi_j*triangelarea;%./3;
 A(j,i)=A(i,j);
 A(i,k)=A(i,k)+gradphi_i'*gradphi_k*triangelarea;%./3;
 A(k,i)=A(i,k);
 A(j,j)=A(j,j)+gradphi_j'*gradphi_j*triangelarea;%./3;
 A(k,k)=A(k,k)+gradphi_k'*gradphi_k*triangelarea;%./3;
 A(j,k)=A(j,k)+gradphi_j'*gradphi_k*triangelarea;%./3;
 A(k,j)=A(j,k);
 A=sparse(A);
 

b(i) = b(i) + fnodi.*triangelarea./3;
b(j) = b(j) + fnodj.*triangelarea./3;
b(k) = b(k) + fnodk.*triangelarea./3;

  
  %A=sparse(A);
end
T=sum(Triangelarea)
A=sparse(A);
b=b;