Home > Src > Level_1 > normah1_1d.m

normah1_1d

PURPOSE ^

NORMAH1_1D Computes H1-norm in 1D

SYNOPSIS ^

function [err_h1]=normah1_1d(fdq, nq, errtype, u, uex, uexx,x, wx, dx, xx, jacx, xy,nov)

DESCRIPTION ^

 NORMAH1_1D   Computes H1-norm in 1D

  [err_h1]=normah1_1d(fdq, nq, errtype, u, uex, uexx,...
            x, wx, dx, xx, jacx, xy,nov)

 Input : fdq = 0  uses Legendre Gauss quadrature formulas with nq nodes
               in each element (exactness degree = 2*nq+1)
             = 1   uses Legendre Gauss Lobatto quadrature formulas on the 
               nodes in which the discrete function is knonwn.
         nq = nodes (in each elemetn) for GL quadrature formulas. Not used if
              fdq=1
         errtype = 0 for absolute error ||u-u_ex||
                   1 for relative error ||u-u_ex||/||u_ex||
         u = numerical solution, column array of length noe
         uex  = exact solution (uex=@(x)[uex(x)], with .*, .^, ./)
         uexx  = first derivative of exact solution 
                 (uexx=@(x)[uex(xx)], with .*, .^, ./)
         x = column array  with LGL nodes in [-1,1]
         wx= column array  with LGL weigths in [-1,1]
         dx= derivative matrix produced with derlgl
         xx = 2-indexes array of size (2,ne): xx(1:2,ie)=[xa_ie;xb_ie]
         jacx = column array of length = ne, containing 
                jacobians of the maps F_ie:[-1,1]---->[xa_ie,xb_ie]
         xy = mesh, column array of length noe
         nov = local to global map. 2-indexes array, size(nov)=[nov,ne]

 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:

SUBFUNCTIONS ^

SOURCE CODE ^

0001 function [err_h1]=normah1_1d(fdq, nq, errtype, u, uex, uexx,...
0002  x, wx, dx, xx, jacx, xy,nov)
0003 % NORMAH1_1D   Computes H1-norm in 1D
0004 %
0005 %  [err_h1]=normah1_1d(fdq, nq, errtype, u, uex, uexx,...
0006 %            x, wx, dx, xx, jacx, xy,nov)
0007 %
0008 % Input : fdq = 0  uses Legendre Gauss quadrature formulas with nq nodes
0009 %               in each element (exactness degree = 2*nq+1)
0010 %             = 1   uses Legendre Gauss Lobatto quadrature formulas on the
0011 %               nodes in which the discrete function is knonwn.
0012 %         nq = nodes (in each elemetn) for GL quadrature formulas. Not used if
0013 %              fdq=1
0014 %         errtype = 0 for absolute error ||u-u_ex||
0015 %                   1 for relative error ||u-u_ex||/||u_ex||
0016 %         u = numerical solution, column array of length noe
0017 %         uex  = exact solution (uex=@(x)[uex(x)], with .*, .^, ./)
0018 %         uexx  = first derivative of exact solution
0019 %                 (uexx=@(x)[uex(xx)], with .*, .^, ./)
0020 %         x = column array  with LGL nodes in [-1,1]
0021 %         wx= column array  with LGL weigths in [-1,1]
0022 %         dx= derivative matrix produced with derlgl
0023 %         xx = 2-indexes array of size (2,ne): xx(1:2,ie)=[xa_ie;xb_ie]
0024 %         jacx = column array of length = ne, containing
0025 %                jacobians of the maps F_ie:[-1,1]---->[xa_ie,xb_ie]
0026 %         xy = mesh, column array of length noe
0027 %         nov = local to global map. 2-indexes array, size(nov)=[nov,ne]
0028 %
0029 % Reference: CHQZ2 = C. Canuto, M.Y. Hussaini, A. Quarteroni, T.A. Zang,
0030 %                    "Spectral Methods. Fundamentals in Single Domains"
0031 %                    Springer Verlag, Berlin Heidelberg New York, 2006.
0032 
0033 %   Written by Paola Gervasio
0034 %   $Date: 2007/04/01$
0035 
0036 
0037 npdx=length(x); [ldnov,ne]=size(nov); noe=nov(npdx,ne);
0038 num=0; den=0;
0039 err_h1=0;
0040 
0041 if fdq==0
0042 
0043 % Legendre Gauss nodes and weigths
0044 
0045 [xg,wxg] = xwlg(nq,-1,1);
0046 
0047 % Evaluation of Lagrange basis polynomials at quadrature nodes xg
0048 %
0049 [phix]= intlag_lgl (x,xg);
0050 
0051 % Loop on spectral elements
0052 
0053 for ie=1:ne
0054 u_loc=u(nov(1:npdx,ie));
0055 [norma_loc1,norma_loc2]=normah1_ex_loc(errtype,npdx,uex,uexx,...
0056  u_loc,x,dx,jacx(ie),xx(1:2,ie),xg,wxg,nq,phix);
0057 num=num+norma_loc1;
0058 den=den+norma_loc2;
0059 end
0060 
0061 elseif fdq==1
0062     
0063 for ie=1:ne
0064 u_loc=u(nov(1:npdx,ie));
0065 xy_loc=xy(nov(1:npdx,ie));
0066 [norma_loc1,norma_loc2]=normah1_loc(errtype,npdx,uex,uexx,...
0067  u_loc,xy_loc, wx,dx,jacx(ie));
0068 num=num+norma_loc1;
0069 den=den+norma_loc2;
0070 end
0071 
0072 end
0073 if errtype==0
0074     err_h1=sqrt(num);
0075 elseif errtype==1
0076     if abs(den)>1.d-14; err_h1=sqrt(num/den); end
0077 end
0078 
0079 return
0080 
0081 function [norma_loc1,norma_loc2]=normah1_ex_loc(errtype,npdx,uex,uexx,...
0082  u, x,dx,jacx,xx,xg,wxg,nq,phix);
0083 
0084 % High degree Legendre Gaussian  formulas to compute H1-norm error
0085 
0086 % mapping quadrature nodes on element ie
0087 
0088 xxg=xg*jacx+(xx(2)+xx(1))*.5;
0089 
0090 % evaluation of exact solution at quadrature nodes.
0091 
0092 U=uex(xxg);
0093 UX=uexx(xxg);
0094 
0095 % Compute numerical derivative
0096 ux=dx*u/jacx; 
0097 
0098 % evaluate numerical solution and its derivative at quadrature nodes.
0099 
0100 u_i=phix*u;
0101 ux_i=phix*ux;
0102 %  figure(1)
0103 %  mesh(xxg,yyg,reshape(U,npx,npy));
0104 %  figure(2)
0105 %  mesh(xxg,yyg,uloc_i)
0106 
0107 % compute the sum
0108 
0109 norma_loc1=sum(((U-u_i).^2+(UX-ux_i).^2).*wxg)*jacx;
0110 
0111 if errtype==0
0112     norma_loc2=0;
0113 else
0114     norma_loc2=sum((U.^2+UX.^2).*wxg)*jacx;
0115 end
0116 
0117 return
0118 
0119 function [norma_loc1,norma_loc2]=normah1_loc(errtype,npdx,uex,uexx, u_loc,...
0120  xy_loc, wx,dx, jacx); 
0121 
0122 % LGL quadrature formulas on npdx nodes to compute H1-norm error
0123 
0124 
0125 % evaluation of exact solution at quadrature nodes.
0126 
0127 U=uex(xy_loc);
0128 UX=uexx(xy_loc);
0129 
0130 % Compute numerical derivative
0131 ux=dx*u_loc/jacx;
0132 
0133 % compute the sum
0134 
0135 norma_loc1=sum(((U-u_loc).^2+(UX-ux).^2).*wx)*jacx;
0136 if errtype==0
0137     norma_loc2=0;
0138 else
0139     norma_loc2=sum((U.^2+UX.^2).*wx)*jacx;
0140 end
0141 
0142 return
0143

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