call_lap_3d

 CALL_LAP_3D Script for pre- and post-processing lap_3d

0001 % CALL_LAP_3D Script for pre- and post-processing lap_3d
0002 %
0003 
0004 %   Written by Paola Gervasio
0005 %   $Date: 2007/04/01$
0006 
0007 %
0008 [uex,uex_x,uex_y,uex_z,ff,gam]=setfun_lap_3d;  % -Delta u + gam *u =f
0009 %
0010 xa=-1;xb=1;   % Omega=(xa,xb) x (ya,yb) x (za,zb)
0011 ya=-1;yb=1;
0012 za=-1;zb=1;
0013 errinf=[];errh1=[];errl2=[];gammax=[]; gammay=[]; gammaz=[];
0014 H=[];
0015 for nex=4;
0016     ney=nex;  % decomposition of Omega in nex * ney * nez parallelepipeds
0017     nez=nex; 
0018     H=[H;2/nex];
0019 for nx=2:4      % polynomial degree in each element along x-direction
0020     ny=nx;    % polynomial degree in each element along y-direction
0021     nz=nx;    % polynomial degree in each element along y-direction
0022     param=zeros(20,1);  
0023     param(1)=1;    % 1=SEM-NI,
0024     param(2)=0;    % 0=no reordering, 1=CM ordering, 2=AMD ordering
0025     param(3)=2;    % 1= solve linear system by \
0026                    % 2= compute extrema eigenvalues of A
0027                    % 3 solve by Schur complement
0028                    % 4= compute extrema eigenvalues of the Schur complement
0029                    % 5= plot eigenvalues of matrix A
0030                    % 6: Compute extrema eigenvalues of  M^{-1}A
0031     param(4)=1;    % 1 computes errors, 0 no
0032     param(5)=0;    % 0 = errors with exact norms, 1= errors with discrete norms
0033     param(6)=nx*2; % n. of nodes for LG quadrature formulas
0034     param(7)=0;    % 0 =absolute errors, 1=relative errors
0035     param(8)=0;    % 0 no plot, 1 mesh, 2 surf, 3 contour
0036     param(9)=(nx+1); % nodes used to plot numerical solution
0037 %     gammax=[0.2]; gammay=[0.2]; % if SEM decomposition is not uniform:
0038 %     gammaz=[0.2];            % they are the arrays with intefaces positions
0039 
0040 % call lap_3d
0041 
0042 [xy,un,D,param]=lap_3d(xa,xb,ya,yb,za,zb,gam,...
0043           uex,uex_x,uex_y,uex_z,ff,nex,nx,ney,ny,nez,nz,gammax,gammay,gammaz,...
0044           param);
0045 
0046 % post
0047 
0048 if (param(3)==1 | param(3)==3)
0049 fprintf('nx=%d,nex=%d,err_inf=%11.4e, err_h1=%11.4e,err_l2=%11.4e \n',...
0050     nx,nex,param(29),param(30),param(31))
0051 errinf=[errinf;param(29)];
0052 errh1=[errh1;param(30)];
0053 errl2=[errl2;param(31)];
0054 
0055 elseif(param(3)==2)
0056 fprintf('nx=%d,nex=%d,lam_max(A)=%11.4e, lam_min(A)=%11.4e, k(A)=%11.4e \n',...
0057     nx,nex,param(23),param(24),abs(param(23))/abs(param(24)));
0058 elseif(param(3)==4)
0059 fprintf('nx=%d,nex=%d,lam_max(S)=%11.4e, lam_min(S)=%11.4e, k(S)=%11.4e \n',...
0060     nx,nex,param(25),param(26),abs(param(25))/abs(param(26)));
0061 elseif(param(3)==6)
0062 fprintf('nx=%d,nex=%d,lam_max(M^{-1}A)=%11.4e, lam_min(M^{-1}A)=%11.4e, k(M^{-1}A)=%11.4e \n',...
0063     nx,nex,param(27),param(28),abs(param(27))/abs(param(28)));
0064 end
0065 
0066 end
0067 end