Home > Src > Level_3 > errors_3d.m

errors_3d

PURPOSE ^

ERRORS_3D Computes errors for 3D b.v.p.

SYNOPSIS ^

function [err_inf,err_h1,err_l2]=errors_3d(x,wx,dx,xx,jacx,y,wy,dy,yy,jacy,z,wz,dz,zz,jacz,xyz,ww,nov,un,uex,uex_x,uex_y,uex_z,param);

DESCRIPTION ^

 ERRORS_3D Computes errors for 3D b.v.p.

 Input:
        x = npdx LGL nodes in [-1,1], previously generated by xwlgl
        wx = npdx LGL weigths in [-1,1], previously generated by xwlgl
        dx = LGL 1-st derivative matrix previously generated by derlgl
        xx = 2-indexes array of size (8,ne)
            xx(1:8,ie)=[x_V1_ie;x_V2_ie;x_V3_ie;x_V4_ie;
                        x_V5_ie;x_V6_ie;x_V7_ie;x_V8_ie]
            (ne=nex*ney is the global number of spectral elements)
        jacx = array (length(jacx)=ne); jacx(ie)= (x_V2_ie-x_V1_ie)/2
        y = npdx LGL nodes in [-1,1], previously generated by xwlgl
        wy = npdx LGL weigths in [-1,1], previously generated by xwlgl
        dy = LGL 1-st derivative matrix previously generated by derlgl
        yy = 2-indexes array of size (8,ne)
            yy(1:8,ie)=[y_V1_ie;y_V2_ie;y_V3_ie;y_V4_ie;
                        y_V5_ie;y_V6_ie;y_V7_ie;y_V8_ie]
            (ne=nex*ney is the global number of spectral elements)
        jacy = array (length(jacy)=ne); jacy(ie)= (y_V3_ie-y_V1_ie)/2
        z = npdz LGL nodes in [-1,1], previously generated by xwlgl
        wz = npdz LGL weigths in [-1,1], previously generated by xwlgl
        dz = LGL 1-st derivative matrix previously generated by derlgl
        zz = 2-indexes array of size (8,ne)
            zz(1:8,ie)=[z_V1_ie;z_V2_ie;z_V3_ie;z_V4_ie;
                        z_V5_ie;z_V6_ie;z_V7_ie;z_V8_ie]
            (ne=nex*ney is the global number of spectral elements)
        jacz = array (length(jacz)=ne); jacz(ie)= (z_V5_ie-z_V1_ie)/2
        xyz = column array with global mesh, length: noe=nov(npdx*npdy*npdz,ne)
        ww = column array with global weigths, 
               length: noe=nov(npdx*npdy*npdz,ne)
              diag(ww) is the mass matrix associated to SEM discretization
        nov = local -global map, previously generated by mesh_3d
        un = numerical solution
        uex  = exact solution (uex=@(x,y,z)[uex(x,y,z)], with .*, .^, ./)
        uex_x  = exact first x-derivative (uex_x=@(x,y,z)[uex_x(x,y,z)], with .*, .^, ./)
        uex_y  = exact first y-derivative (uex_y=@(x,y,z)[uex_y(x,y,z)], with .*, .^, ./)
        uex_z  = exact first z-derivative (uex_z=@(x,y,z)[uex_z(x,y,z)], with .*, .^, ./)
        param(4) = 1: compute errors (L^inf-norm, L2-norm, H1-norm)
                      on the exact solution
                   2: no  errors are computed
        param(5) = 0: LG quadrature formulas with high precision degree are
                      used to compute norms (exact norms)
                   1: LGL quadrature formulas with npdx,npdy,npdz nodes are
                      used to compute norms (discrete norms)
                   (used only if param(4) == 1)
        param(6) = number of nodes for high degree quadrature formula,
                   (used only if param(5) == 0 & param(4) == 1)
        param(7) = 0: absolute errors are computed
                   1: relative errors are computed
                   (used only if param(4) == 1)

 Output:

         err_inf = ||u_ex-u_N|| with respect to discrete maximum-norm
         err_h1 = ||u_ex-u_N|| with respect to discrete H1-norm
         err_l2 = ||u_ex-u_N|| with respect to discrete L2-norm

 Reference: CHQZ2 = C. Canuto, M.Y. Hussaini, A. Quarteroni, T.A. Zang,
                    "Spectral Methods. Fundamentals in Single Domains"
                    Springer Verlag, Berlin Heidelberg New York, 2006.

CROSS-REFERENCE INFORMATION ^

This function calls: This function is called by:

SOURCE CODE ^

0001 function [err_inf,err_h1,err_l2]=errors_3d(x,wx,dx,xx,jacx,y,wy,dy,yy,jacy,...
0002 z,wz,dz,zz,jacz,...
0003 xyz,ww,nov,un,uex,uex_x,uex_y,uex_z,param);
0004 % ERRORS_3D Computes errors for 3D b.v.p.
0005 %
0006 % Input:
0007 %        x = npdx LGL nodes in [-1,1], previously generated by xwlgl
0008 %        wx = npdx LGL weigths in [-1,1], previously generated by xwlgl
0009 %        dx = LGL 1-st derivative matrix previously generated by derlgl
0010 %        xx = 2-indexes array of size (8,ne)
0011 %            xx(1:8,ie)=[x_V1_ie;x_V2_ie;x_V3_ie;x_V4_ie;
0012 %                        x_V5_ie;x_V6_ie;x_V7_ie;x_V8_ie]
0013 %            (ne=nex*ney is the global number of spectral elements)
0014 %        jacx = array (length(jacx)=ne); jacx(ie)= (x_V2_ie-x_V1_ie)/2
0015 %        y = npdx LGL nodes in [-1,1], previously generated by xwlgl
0016 %        wy = npdx LGL weigths in [-1,1], previously generated by xwlgl
0017 %        dy = LGL 1-st derivative matrix previously generated by derlgl
0018 %        yy = 2-indexes array of size (8,ne)
0019 %            yy(1:8,ie)=[y_V1_ie;y_V2_ie;y_V3_ie;y_V4_ie;
0020 %                        y_V5_ie;y_V6_ie;y_V7_ie;y_V8_ie]
0021 %            (ne=nex*ney is the global number of spectral elements)
0022 %        jacy = array (length(jacy)=ne); jacy(ie)= (y_V3_ie-y_V1_ie)/2
0023 %        z = npdz LGL nodes in [-1,1], previously generated by xwlgl
0024 %        wz = npdz LGL weigths in [-1,1], previously generated by xwlgl
0025 %        dz = LGL 1-st derivative matrix previously generated by derlgl
0026 %        zz = 2-indexes array of size (8,ne)
0027 %            zz(1:8,ie)=[z_V1_ie;z_V2_ie;z_V3_ie;z_V4_ie;
0028 %                        z_V5_ie;z_V6_ie;z_V7_ie;z_V8_ie]
0029 %            (ne=nex*ney is the global number of spectral elements)
0030 %        jacz = array (length(jacz)=ne); jacz(ie)= (z_V5_ie-z_V1_ie)/2
0031 %        xyz = column array with global mesh, length: noe=nov(npdx*npdy*npdz,ne)
0032 %        ww = column array with global weigths,
0033 %               length: noe=nov(npdx*npdy*npdz,ne)
0034 %              diag(ww) is the mass matrix associated to SEM discretization
0035 %        nov = local -global map, previously generated by mesh_3d
0036 %        un = numerical solution
0037 %        uex  = exact solution (uex=@(x,y,z)[uex(x,y,z)], with .*, .^, ./)
0038 %        uex_x  = exact first x-derivative (uex_x=@(x,y,z)[uex_x(x,y,z)], with .*, .^, ./)
0039 %        uex_y  = exact first y-derivative (uex_y=@(x,y,z)[uex_y(x,y,z)], with .*, .^, ./)
0040 %        uex_z  = exact first z-derivative (uex_z=@(x,y,z)[uex_z(x,y,z)], with .*, .^, ./)
0041 %        param(4) = 1: compute errors (L^inf-norm, L2-norm, H1-norm)
0042 %                      on the exact solution
0043 %                   2: no  errors are computed
0044 %        param(5) = 0: LG quadrature formulas with high precision degree are
0045 %                      used to compute norms (exact norms)
0046 %                   1: LGL quadrature formulas with npdx,npdy,npdz nodes are
0047 %                      used to compute norms (discrete norms)
0048 %                   (used only if param(4) == 1)
0049 %        param(6) = number of nodes for high degree quadrature formula,
0050 %                   (used only if param(5) == 0 & param(4) == 1)
0051 %        param(7) = 0: absolute errors are computed
0052 %                   1: relative errors are computed
0053 %                   (used only if param(4) == 1)
0054 %
0055 % Output:
0056 %
0057 %         err_inf = ||u_ex-u_N|| with respect to discrete maximum-norm
0058 %         err_h1 = ||u_ex-u_N|| with respect to discrete H1-norm
0059 %         err_l2 = ||u_ex-u_N|| with respect to discrete L2-norm
0060 %
0061 % Reference: CHQZ2 = C. Canuto, M.Y. Hussaini, A. Quarteroni, T.A. Zang,
0062 %                    "Spectral Methods. Fundamentals in Single Domains"
0063 %                    Springer Verlag, Berlin Heidelberg New York, 2006.
0064 
0065 %   Written by Paola Gervasio
0066 %   $Date: 2007/04/01$
0067 
0068 
0069 fdq=param(5); 
0070 if fdq ==0
0071 nq=param(6);
0072 else 
0073 nq=length(wx);
0074 end
0075 errtype=param(7);
0076 
0077 % compute errors on the exact solution
0078 
0079 npdx=length(wx); npdy=length(wy); npdz=length(wz);
0080 [ldnov,ne]=size(nov);
0081 noe=length(un);
0082 
0083 % Evaluate exact solution
0084 
0085 u=zeros(noe,1)+uex(xyz(:,1),xyz(:,2),xyz(:,3));
0086 
0087 %  compute difference
0088 
0089 err=u-un;
0090 
0091 % ||u_ex-u_N|| with respect to discrete maximum-norm
0092 
0093 err_inf=norm(err,inf);
0094 if errtype~=0
0095     err_inf=err_inf/norm(u,inf);
0096 end
0097 % ||u_ex-u_N|| with respect to H1 norm
0098 
0099 % fdq=0 LG quadrature formula with nq nodes on each element
0100 % fdq=1 LGL quadrature formula with npdx nodes on each element
0101 
0102 [err_h1]=normah1_3d(fdq,nq,errtype,x,wx,dx,xx,jacx,y,wy,dy,yy,jacy,...
0103 z,wz,dz,zz,jacz,xyz,ww,nov,un,u,uex,uex_x,uex_y,uex_z);
0104 
0105 % ||u_ex-u_N|| with respect to L2 norm
0106 
0107 [err_l2]=normal2_3d(fdq,nq,errtype,x,wx,xx,jacx,y,wy,yy,jacy,...
0108 z,wz,dz,zz,jacz,xyz,ww,nov,un,u,uex);
0109 
0110 return

Generated on Fri 21-Sep-2007 10:07:00 by m2html © 2003