薛定谔特征值问题和拉普拉斯特征值问题的有限元MATLAB程序（一维）
−(p(x)u′(x))′+q(x)u(x)=λω(x)u(x)
in Ω=(a,b),u(a)=0,u(b)=0
这是对一般的Sturm-Liouville问题的一个表述。下面主要谈这个问题的两个特例，拉普拉斯特征值问题和薛定谔特征值问题。
拉普拉斯特征值问题是：
−u′′(x)=λu(x) in Ω=(0,1)u(0)=0,u(1)=0
作为SL问题的一个特例，它有精确解：
un(x)=sin(nπx),λn=n2π2(n=1,2...)带调和势的薛定谔特征值问题的一般表述是：
−u′′(x)2+V(x)u(x)=ϵu(x) in Ω=(−L,−L)u(−L)=0,u(L)=0
这里的L在程序中给定为10， V(x) 取为 x22 。略过简单的变分过程，下面直接给出这两类问题的一个程序：
两个调用函数的主程序：
拉普拉斯主程序：
clear all; close all;
options = struct('verbose','on','problem','laplace','domain',[0 1], ...
'analysis','fe', ...
'element_type','spectral','element_order',5,'numel',8, ...
'nsp',7, ... % number of quadrature points
'eigsolver','eigs','eigtol',1e-14,'nbeigs',6);
[lambda,numdofs,cond] = fem1d_evp(options);
sumeig = sum(lambda);
% Exact solution = n^2 pi^2 , n = 1, 2, ...
exact = 91*pi^2; % 1 + 4 + 9 + 16 + 25 + 36 = 91
erreig = sumeig - exact;
fprintf('Error in the sum of the first six eigenvalues = %15.7e\n',erreig);
薛定谔主程序：
clear all; close all;
options = struct('verbose','on','problem','schrodinger','domain',[-10 10], ...
'analysis','fe', ...
'element_type','spectral','element_order',5,'numel',8, ...
'nsp',7, ... % number of quadrature points
'eigsolver','eigs','eigtol',1e-14,'nbeigs',6);
[lambda,numdofs,cond] = fem1d_evp(options);
sumeig = sum(lambda);
% Exact solution = n + 1/2 , n = 0,1,....
exact = 18; % 1/2 + 3/2 + 5/2 + 7/2 + 9/2 + 11/2 = 18
erreig = sumeig - exact;
fprintf('Error in the sum of the first six eigenvalues = %15.7e\n',erreig);
下面给出几个主要的子程序：
- fem1d_evp.m求解特征值问题的主要函数：
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% Scientific Computing : Computational Methods for Nonlinear
% Eigenvalue Problems, LSEC, Chinese Academy of Sciences, Beijing, China.
%============================================================================
%
% Finite Element Code for Laplace and Schrodinger Eigenvalue Problems
% -------------------------------------------------------------------
% Author : Song Lu, LSEC
% E-mail : lusongcool@163.com
%
%============================================================================
%
% MODEL 1D Problem (Strong Form)
% ==============================
% One-dimensional spectral finite element program to solve the Laplace and
% Schrodinger eigenproblems:
%
% -u''/A + V u = \lambda u in \Omega = (a,b),
%
% where V is the harmonic potential in the Schrodinger equation and V = 0 for
% the Laplace eigenproblem). The constant A is 2 for the Schrodinger and 1
% for the Laplace problem. \lambda is the eigenvalue and u is the eigenfunction
%
% Boundary conditions are: u(0) = u(1) = 0 for Laplace eigenproblem
% u(-10) = u(10) = 0 for Schrodinger eigenproblem
%
% FE Approximation
% ================
% Galerkin method with spectral finite elements
%
% -------------------------------------------------------
% USAGE: fem1d_evp(options)
% -------------------------------------------------------
% options : structure
% : options.verbose = 'on' or 'off';
% : options.problem = 'laplace' or 'schrodinger'
% : options.domain = [0 1] or [-10 10], etc.
% : options.numel = number of elements
% : options.nsp = number of quadrature points
% : options.analysis = 'fe'
% Default setting :
% options = struct('verbose','on','problem','laplace','domain',[0 1], ...
% 'analysis','fe','element_type','spectral', ...
% 'element_order',4,'numel',3,'nsp',6, ...
% 'eigsolver','eigs','eigtol',1e-12,'nbeigs',6);
%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
function [lambda,numdofs,condbigs] = fem1d_evp(options)
%Hu = lambdaSu
global bigh;%H
global bigs;%S
condbigs = 0;%H的条件数
if nargin < 1
help(mfilename);
options = struct('verbose','on','problem','laplace','domain',[0 1],'analysis','fe',...
'element_type','spectral','element_order',4,'numel',3,...
'nsp',6,'eigsolver','eigs','eigtol',1e-12,'nbeigs',6);
end
addpath(genpath(pwd));%将子文件的路径
problem = options.problem;
if (~strcmp(problem,'laplace') && ~strcmp(problem,'schrodinger'))
error('Aborting: Problem name not known');
end%主要考虑拉普拉斯方程和薛定谔方程
numel = options.numel;
order = options.element_order;%每个基元分几段
nsp = options.nsp;
numnod = order*numel + 1;%节点总数
options.numnod = numnod;
%
% start timing
%
tStart = tic;%计时函数
%
% call mesher
%
fprintf('Constructing mesh . . .\n\n');
[coord,connect] = mesher(options);%构建洛巴托积分点坐标和连接矩阵
%coord(2:end,:)-coord(1:end-1,:) coordy洛巴托积分点坐标
fprintf('Arrays for dof info \n\n');
[nodaldofs,dofinfo,numdofs] = setupdofs(options,coord,connect);
%找出哪些分点是自由的，nodaldofs是结构体，指出每个点是第几个自由点，dofinfo是
%自由点在整个节点中的编号，numdofs是自由点编号。
options.numdofs = numdofs;
fprintf('Allocate arrays for I, J, XH, XS: sparse assembly\n\n');
m = 0; % size of index arrays
for e = 1:numel
conn = connect(e,:); sctr = [ nodaldofs(conn).list ];
m = m + (length(sctr))^2;
end
%根据自由度，计算来存放各个矩阵需要的空间
I = zeros(m,1); J = I; XH = I; XS = I;%初始化
%
% compute and assemble Hamiltonian matrix H and overlap matrix S
%
% loop over all the elements
%
ind = 1; % indices for arrays I, J, XH, XS
for e = 1:numel
printelementnumber(e,numel);%根据基元数情况打印编号
conn = connect(e,:); sctr = [ nodaldofs(conn).list ];
%提取每个单元上节点的整体序号，以及自由点的内部序数
%
% compute element contributions to Hamiltonian H (kinectic energy) and overlap S
%
[he,se] = elemhamiltonianoverlap(options,e,coord,conn,sctr,dofinfo);
%组装单元刚度矩阵
% insert he, se into vector arrays
%
sctrlen = length(sctr);
temp = meshgrid(sctr);
I(ind:ind+sctrlen^2-1) = temp(:);%把temp按列拉直
temp = temp';
J(ind:ind+sctrlen^2-1) = temp(:);%按行拉直
XH(ind:ind+sctrlen^2-1) = he(:);
XS(ind:ind+sctrlen^2-1) = se(:);
ind = ind + sctrlen^2;
end
%
% Create H and S-matrices from the arrays
%
fprintf('Performing sparse assembly of system matrices (%d x %d)\n\n',numdofs,numdofs);
bigh = sparse(I,J,XH,numdofs,numdofs);
bigs = sparse(I,J,XS,numdofs,numdofs);
fprintf('Freeing memory: I, J, XH, XS arrays\n\n',numdofs,numdofs);
clear I; clear J; clear XH; clear XS;
condbigs = condest(bigs);
fprintf('Condition number of S = %g\n',condbigs);
%
% solve the generalized eigenvalue problem: H c = \lambda S c
%
fprintf('Solving generalized eigenproblem . . .\n\n');
[lambda,eigv] = eigsolverinterface(options.eigsolver,options.nbeigs);
%%%%%%%%%%%%%%%%%%
% DISPLAY RESULTS
%%%%%%%%%%%%%%%%%%
fprintf('\n\nFEM RESULTS\n');
fprintf('===========\n');
format long;
fprintf('Problem domain = [%g %g]\n', options.domain(1), options.domain(2));
fprintf('Number of elements = %d\n', options.numel);
fprintf('Number of nodes = %d\n', options.numnod);
fprintf('Number of quadrature points = %d\n', nsp);
fprintf('Number of DOFs = %d\n', options.numdofs);
fprintf('Condition number of S = %g\n',condbigs);
fprintf('Eigenvalues = ');
fprintf('%18.10f ', sort(lambda)); % ascending order
fprintf('\n\n');
tEnd = toc(tStart);
fprintf('Total time in fem1d_evp: %d minutes and %f seconds\n\n',floor(tEnd/60),rem(tEnd,60));
return
end
- mesher.m定义节点以及各点间的相互关系（连接矩阵）
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% Purpose
% =======
% Finite element mesh in 1D with option to create regular (m-node) elements.
% Output plot is also generated to show the mesh.
%
%----------------------------------------------------------------------------
% Function Syntax
% ===============
% function [coord,connect] = mesher(options)
%
% Inputs
% ======
% options : data structure with mesh/element information
%
% Outputs
% =======
% coord : 1D array which stores the coordinates of the nodes
% coord(I,1) : node I coordinates
% connect : element connectivity
%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%
function [coord,connect] = mesher(options)
element_type = options.element_type;
order = options.element_order;
switch element_type
case 'spectral'
otherwise
fprintf('Mesh generator not coded for element type = %s'\n,element_type);
error('Aborting . .');
end
numel = options.numel;
plot = options.verbose;
domain = options.domain;
xmin = domain(1); xmax = domain(2);
numnod = order*numel+1;%节点数目
coord = zeros(numnod,1);
coord(:,1) = linspace(xmin,xmax,numnod);%坐标等分
%
% connectivity
%
connect = zeros(numel,order+1);%连接矩阵，没行表示一个基元块
for e = 1:numel
is = order*(e-1) + 1;
connect(e,:) = is:is+order;
end
%
% modify nodal coordinates (Lobatto points)
%
if (order >= 2)
xi = lobatto_points(order+1);% [-1 1]上的洛巴托积分点
for e = 1:numel
is = order*(e-1) + 1;
a = coord(is,:); b = coord(is+order,:);%不同基元上的左右两点
for n = is:is+order
index = n-is+1;%节点的元内编号
coord(n,:) = 0.5*(b + a) + 0.5*(b - a).*xi(index);
%通过映射关系找出每个基元上的洛巴托点对应的点
end
end
end
%把均分点调整成洛巴托积分点
%-------------------------------------------------------------------
% Mesh plotting
%-------------------------------------------------------------------
if (strcmp(plot,'on'))
hold on;
axis([xmin xmax -1 1 ]);
xlabel('x','fontsize', 12); % label for x-axis
plotnodes(coord);
hold off;
end
return
end
- plotnodes.m画一下修正为洛巴托积分点的各个节点，其实不是必要的：
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% Purpose
% =======
% Plot nodes
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
function plotnodes(coord)
close all;
x = coord(:,1);
y = zeros(1,length(x));
set(gcf,'Position',[0 0 600 100]);
plot(x,y,'.','Linewidth',3,'MarkerEdgeColor','b','MarkerFaceColor','b','MarkerSize',25);
return
end
- eigsolverinterface.m用于求解矩阵特征值问题的特征值和特征向量，比起eigs函数，多了一步对非对称矩阵做对称化操作，大大加快了求解速度。
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
function [lambda,eigv] = eigsolverinterface(solvertype,nbeigs)
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
global bigh;
global bigs;
solvertype = lower(solvertype);
%统统变小写
%
% solve the generalized eigenvalue problem: H c = \lambda S c
%
nbfreedofs = length(bigh(:,1));
switch solvertype
case 'eigs'
opts.maxit = 1000;
opts.tol = 1.e-12;
opts.disp = 2;
% symmetrize
bigs = 0.5*(bigs + bigs');
[eigv,D,failureflag] = eigs(bigh,bigs,nbeigs,'SM',opts);
for i = 1:nbeigs
lambda(i,1) = D(i,i);
end
otherwise
error('Aborting - solvertype undefined');
end
- elemhamiltonianoverlap.m计算每个单元上的刚度和载荷矩阵。
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% function [he,se] = elemhamiltonianoverlap(options,e,coord,nodes,dofs,dofinfo)
% Purpose
% =======
% Element Hamiltonian and overlap matrices for spectral finite elements
%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
function [he,se] = elemhamiltonianoverlap(options,e,coord,nodes,dofs,dofinfo)
% initialize element Hamiltonian and overlap matrices to zero
%
nsp = options.nsp;%高斯点数
problem = options.problem;%
numel = options.numel;%基元数目
fac = 0.5;%针对于薛定谔问题的系数
if (strcmp(problem,'laplace')) fac = 1.0; end
eorder = options.element_order;%order就是每个单元涉及到的基函数的数目
dofs_to_nodes = dofinfo(dofs); % dofs = DOFs of the element，当前单元自由节点对应的总节点编号
ndofs = length(dofs_to_nodes); % global node Ids that correspond to these DOFs，这个单元自由节点的数目
dofs_to_locindices = zeros(1,ndofs);
for i = 1:ndofs
dofs_to_locindices(i) = find(nodes == dofs_to_nodes(i));
%自由节点们在单元中的序号
end
xe = coord(nodes,1); % xe是单元节点的洛巴托坐标
a = xe(1); b = xe(end);%单元两端坐标（洛巴托）
Le = b - a; detj = Le/2.;%区间长度与伸缩比
he = zeros(ndofs,ndofs);
se = zeros(ndofs,ndofs);%单刚初始化
%
% loop over Gauss points in xi-direction
%
gauss = gauss_points(nsp); % Using Gauss-Legendre quadrature (not Gauss-Lobatto !)
weight = gauss_weights(nsp);%获取6个高斯点以及权重
for i = 1:nsp %按高斯点来对基函数配对
xi = gauss(i);
x = 0.5*(1. - xi)*a + 0.5*(1. + xi)*b;%x为高斯点在该段区间中对应的值
%
% FE shape functions and local derivatives
%
[NFE, dNFEdxi ] = shapefunctions(eorder,xi); % all the shapefunctions
%各个标准区间上的节点基函数在高斯点出的函数值
N = NFE(dofs_to_locindices);
dNdxi = dNFEdxi(dofs_to_locindices);
%高斯点在具有自由度的几个基函数上的值
% derivative of shape function with respect to global coordinate x
dNdx = dNdxi/detj;
%导数值映射到一般区间上会变
%
% assemble element Hamiltonian (he) and overlap matrices (se)
%
he = he + fac*weight(i)*(dNdx'*dNdx)*detj;
se = se + weight(i)*N'*N*detj;
if (strcmp(problem,'schrodinger'))
% V = x^2/2 - Harmonic potential
he = he + 0.5*x*x*weight(i)*N'*N*detj;
end
end
return
end
- gauss_points.m,gauss_weights.m,lobatto_points.m,lobatto_weights.m 顾名思义，分别给出高斯点和洛巴托点及其权重：
function x = guass_points(n);
%
% function x = gauss_points(n)
%
% For 1 <= n <= 30, returns the abscissas x of an n point
% Gauss-Legendre quadrature rule over the interval [-1,1].
%
x = ones(1,n);
if ( n == 1 )
x(1) = 0.0;
elseif ( n == 2 )
x(1) = - 0.577350269189625764509148780502;
x(2) = 0.577350269189625764509148780502;
elseif ( n == 3 )
x(1) = - 0.774596669241483377035853079956;
x(2) = 0.000000000000000;
x(3) = 0.774596669241483377035853079956;
elseif ( n == 4 )
x(1) = - 0.861136311594052575223946488893;
x(2) = - 0.339981043584856264802665759103;
x(3) = 0.339981043584856264802665759103;
x(4) = 0.861136311594052575223946488893;
elseif ( n == 5 )
x(1) = - 0.906179845938663992797626878299;
x(2) = - 0.538469310105683091036314420700;
x(3) = 0.0;
x(4) = 0.538469310105683091036314420700;
x(5) = 0.906179845938663992797626878299;
elseif ( n == 6 )
x(1) = - 0.932469514203152027812301554494;
x(2) = - 0.661209386466264513661399595020;
x(3) = - 0.238619186083196908630501721681;
x(4) = 0.238619186083196908630501721681;
x(5) = 0.661209386466264513661399595020;
x(6) = 0.932469514203152027812301554494;
elseif ( n == 7 )
x(1) = - 0.949107912342758524526189684048;
x(2) = - 0.741531185599394439863864773281;
x(3) = - 0.405845151377397166906606412077;
x(4) = 0.0;
x(5) = 0.405845151377397166906606412077;
x(6) = 0.741531185599394439863864773281;
x(7) = 0.949107912342758524526189684048;
elseif ( n == 8 )
x(1) = - 0.960289856497536231683560868569;
x(2) = - 0.796666477413626739591553936476;
x(3) = - 0.525532409916328985817739049189;
x(4) = - 0.183434642495649804939476142360;
x(5) = 0.183434642495649804939476142360;
x(6) = 0.525532409916328985817739049189;
x(7) = 0.796666477413626739591553936476;
x(8) = 0.960289856497536231683560868569;
elseif ( n == 9 )
x(1) = - 0.968160239507626089835576202904;
x(2) = - 0.836031107326635794299429788070;
x(3) = - 0.613371432700590397308702039341;
x(4) = - 0.324253423403808929038538014643;
x(5) = 0.0;
x(6) = 0.324253423403808929038538014643;
x(7) = 0.613371432700590397308702039341;
x(8) = 0.836031107326635794299429788070;
x(9) = 0.968160239507626089835576202904;
elseif ( n == 10 )
x(1) = - 0.973906528517171720077964012084;
x(2) = - 0.865063366688984510732096688423;
x(3) = - 0.679409568299024406234327365115;
x(4) = - 0.433395394129247290799265943166;
x(5) = - 0.148874338981631210884826001130;
x(6) = 0.148874338981631210884826001130;
x(7) = 0.433395394129247290799265943166;
x(8) = 0.679409568299024406234327365115;
x(9) = 0.865063366688984510732096688423;
x(10) = 0.973906528517171720077964012084;
elseif (n == 11)
x(1)=-0.978228658146056992803938001123;
x(2)=-0.887062599768095299075157769304;
x(3)=-0.730152005574049324093416252031;
x(4)=-0.519096129206811815925725669459;
x(5)=-0.269543155952344972331531985401;
x(6)= 0;
x(7)= 0.269543155952344972331531985401;
x(8)= 0.519096129206811815925725669459;
x(9)= 0.730152005574049324093416252031;
x(10)= 0.887062599768095299075157769304;
x(11)= 0.978228658146056992803938001123;
elseif (n == 12)
x(1)=-0.981560634246719250690549090149;
x(2)=-0.904117256370474856678465866119;
x(3)=-0.769902674194304687036893833213;
x(4)=-0.587317954286617447296702418941;
x(5)=-0.367831498998180193752691536644;
x(6)=-0.125233408511468915472441369464;
x(7)= 0.125233408511468915472441369464;
x(8)= 0.367831498998180193752691536644;
x(9)= 0.587317954286617447296702418941;
x(10)= 0.769902674194304687036893833213;
x(11)= 0.904117256370474856678465866119;
x(12)= 0.981560634246719250690549090149;
elseif (n == 13)
x(1)=-0.984183054718588149472829448807;
x(2)=-0.917598399222977965206547836501;
x(3)=-0.801578090733309912794206489583;
x(4)=-0.642349339440340220643984606996;
x(5)=-0.448492751036446852877912852128;
x(6)=-0.230458315955134794065528121098;
x(7)= 0;
x(8)= 0.230458315955134794065528121098;
x(9)= 0.448492751036446852877912852128;
x(10)= 0.642349339440340220643984606996;
x(11)= 0.801578090733309912794206489583;
x(12)= 0.917598399222977965206547836501;
x(13)= 0.984183054718588149472829448807;
elseif (n == 14)
x(1)=-0.986283808696812338841597266704;
x(2)=-0.928434883663573517336391139378;
x(3)=-0.827201315069764993189794742650;
x(4)=-0.687292904811685470148019803019;
x(5)=-0.515248636358154091965290718551;
x(6)=-0.319112368927889760435671824168;
x(7)=-0.108054948707343662066244650220;
x(8)= 0.108054948707343662066244650220;
x(9)= 0.319112368927889760435671824168;
x(10)= 0.515248636358154091965290718551;
x(11)= 0.687292904811685470148019803019;
x(12)= 0.827201315069764993189794742650;
x(13)= 0.928434883663573517336391139378;
x(14)= 0.986283808696812338841597266704;
elseif (n == 15)
x(1)=-0.987992518020485428489565718587;
x(2)=-0.937273392400705904307758947710;
x(3)=-0.848206583410427216200648320774;
x(4)=-0.724417731360170047416186054614;
x(5)=-0.570972172608538847537226737254;
x(6)=-0.394151347077563369897207370981;
x(7)=-0.201194093997434522300628303395;
x(8)= 0;
x(9)= 0.201194093997434522300628303395;
x(10)= 0.394151347077563369897207370981;
x(11)= 0.570972172608538847537226737254;
x(12)= 0.724417731360170047416186054614;
x(13)= 0.848206583410427216200648320774;
x(14)= 0.937273392400705904307758947710;
x(15)= 0.987992518020485428489565718587;
elseif (n == 16)
x(1)=-0.989400934991649932596154173450;
x(2)=-0.944575023073232576077988415535;
x(3)=-0.865631202387831743880467897712;
x(4)=-0.755404408355003033895101194847;
x(5)=-0.617876244402643748446671764049;
x(6)=-0.458016777657227386342419442984;
x(7)=-0.281603550779258913230460501460;
x(8)=-0.0950125098376374401853193354250;
x(9)= 0.0950125098376374401853193354250;
x(10)= 0.281603550779258913230460501460;
x(11)= 0.458016777657227386342419442984;
x(12)= 0.617876244402643748446671764049;
x(13)= 0.755404408355003033895101194847;
x(14)= 0.865631202387831743880467897712;
x(15)= 0.944575023073232576077988415535;
x(16)= 0.989400934991649932596154173450;
elseif (n == 17)
x(1)=-0.990575475314417335675434019941;
x(2)=-0.950675521768767761222716957896;
x(3)=-0.880239153726985902122955694488;
x(4)=-0.781514003896801406925230055520;
x(5)=-0.657671159216690765850302216643;
x(6)=-0.512690537086476967886246568630;
x(7)=-0.351231763453876315297185517095;
x(8)=-0.178484181495847855850677493654;
x(9)= 0;
x(10)= 0.178484181495847855850677493654;
x(11)= 0.351231763453876315297185517095;
x(12)= 0.512690537086476967886246568630;
x(13)= 0.657671159216690765850302216643;
x(14)= 0.781514003896801406925230055520;
x(15)= 0.880239153726985902122955694488;
x(16)= 0.950675521768767761222716957896;
x(17)= 0.990575475314417335675434019941;
elseif (n == 18)
x(1)=-0.991565168420930946730016004706;
x(2)=-0.955823949571397755181195892930;
x(3)=-0.892602466497555739206060591127;
x(4)=-0.803704958972523115682417455015;
x(5)=-0.691687043060353207874891081289;
x(6)=-0.559770831073947534607871548525;
x(7)=-0.411751161462842646035931793833;
x(8)=-0.251886225691505509588972854878;
x(9)=-0.0847750130417353012422618529358;
x(10)= 0.0847750130417353012422618529358;
x(11)= 0.251886225691505509588972854878;
x(12)= 0.411751161462842646035931793833;
x(13)= 0.559770831073947534607871548525;
x(14)= 0.691687043060353207874891081289;
x(15)= 0.803704958972523115682417455015;
x(16)= 0.892602466497555739206060591127;
x(17)= 0.955823949571397755181195892930;
x(18)= 0.991565168420930946730016004706;
elseif (n == 19)
x(1)=-0.992406843843584403189017670253;
x(2)=-0.960208152134830030852778840688;
x(3)=-0.903155903614817901642660928532;
x(4)=-0.822714656537142824978922486713;
x(5)=-0.720966177335229378617095860824;
x(6)=-0.600545304661681023469638164946;
x(7)=-0.464570741375960945717267148104;
x(8)=-0.316564099963629831990117328850;
x(9)=-0.160358645640225375868096115741;
x(10)= 0;
x(11)= 0.160358645640225375868096115741;
x(12)= 0.316564099963629831990117328850;
x(13)= 0.464570741375960945717267148104;
x(14)= 0.600545304661681023469638164946;
x(15)= 0.720966177335229378617095860824;
x(16)= 0.822714656537142824978922486713;
x(17)= 0.903155903614817901642660928532;
x(18)= 0.960208152134830030852778840688;
x(19)= 0.992406843843584403189017670253;
elseif (n == 20)
x(1)=-0.993128599185094924786122388471;
x(2)=-0.963971927277913791267666131197;
x(3)=-0.912234428251325905867752441203;
x(4)=-0.839116971822218823394529061702;
x(5)=-0.746331906460150792614305070356;
x(6)=-0.636053680726515025452836696226;
x(7)=-0.510867001950827098004364050955;
x(8)=-0.373706088715419560672548177025;
x(9)=-0.227785851141645078080496195369;
x(10)=-0.0765265211334973337546404093988;
x(11)= 0.0765265211334973337546404093988;
x(12)= 0.227785851141645078080496195369;
x(13)= 0.373706088715419560672548177025;
x(14)= 0.510867001950827098004364050955;
x(15)= 0.636053680726515025452836696226;
x(16)= 0.746331906460150792614305070356;
x(17)= 0.839116971822218823394529061702;
x(18)= 0.912234428251325905867752441203;
x(19)= 0.963971927277913791267666131197;
x(20)= 0.993128599185094924786122388471;
elseif (n == 21)
x(1)=-0.99375217062038950026024203593794;
x(2)=-0.96722683856630629431662221490770;
x(3)=-0.92009933415040082879018713371497;
x(4)=-0.85336336458331728364725063858757;
x(5)=-0.76843996347567790861587785130623;
x(6)=-0.66713880419741231930596666999034;
x(7)=-0.55161883588721980705901879672431;
x(8)=-0.42434212020743878357366888854379;
x(9)=-0.28802131680240109660079251606460;
x(10)=-0.14556185416089509093703098233869;
x(11)=0;
x(12)=0.14556185416089509093703098233869;
x(13)=0.28802131680240109660079251606460;
x(14)=0.42434212020743878357366888854379;
x(15)=0.55161883588721980705901879672431;
x(16)=0.66713880419741231930596666999034;
x(17)=0.76843996347567790861587785130623;
x(18)=0.85336336458331728364725063858757;
x(19)=0.92009933415040082879018713371497;
x(20)=0.96722683856630629431662221490770;
x(21)=0.99375217062038950026024203593794;
elseif (n == 22)
x(1)=-0.99429458548239929207303142116130;
x(2)=-0.97006049783542872712395098676527;
x(3)=-0.92695677218717400052069293925905;
x(4)=-0.86581257772030013653642563701938;
x(5)=-0.78781680597920816200427795540835;
x(6)=-0.69448726318668278005068983576226;
x(7)=-0.58764040350691159295887692763865;
x(8)=-0.46935583798675702640633071096641;
x(9)=-0.34193582089208422515814742042738;
x(10)=-0.20786042668822128547884653391955;
x(11)=-0.069739273319722221213841796118628;
x(12)=0.069739273319722221213841796118628;
x(13)=0.20786042668822128547884653391955;
x(14)=0.34193582089208422515814742042738;
x(15)=0.46935583798675702640633071096641;
x(16)=0.58764040350691159295887692763865;
x(17)=0.69448726318668278005068983576226;
x(18)=0.78781680597920816200427795540835;
x(19)=0.86581257772030013653642563701938;
x(20)=0.92695677218717400052069293925905;
x(21)=0.97006049783542872712395098676527;
x(22)=0.99429458548239929207303142116130;
elseif (n == 23)
x(1)=-0.99476933499755212352392571544557;
x(2)=-0.97254247121811523195602407682078;
x(3)=-0.93297108682601610234919698903842;
x(4)=-0.87675235827044166737815688593415;
x(5)=-0.80488840161883989215111840699678;
x(6)=-0.71866136313195019446162448374862;
x(7)=-0.61960987576364615638509731164960;
x(8)=-0.50950147784600754968979304786685;
x(9)=-0.39030103803029083142148887288061;
x(10)=-0.26413568097034493053386953828331;
x(11)=-0.13325682429846611093174268224177;
x(12)=0;
x(13)=0.13325682429846611093174268224177;
x(14)=0.26413568097034493053386953828331;
x(15)=0.39030103803029083142148887288061;
x(16)=0.50950147784600754968979304786685;
x(17)=0.61960987576364615638509731164960;
x(18)=0.71866136313195019446162448374862;
x(19)=0.80488840161883989215111840699678;
x(20)=0.87675235827044166737815688593415;
x(21)=0.93297108682601610234919698903842;
x(22)=0.97254247121811523195602407682078;
x(23)=0.99476933499755212352392571544557;
elseif (n == 24)
x(1)=-0.99518721999702136017999740970074;
x(2)=-0.97472855597130949819839199300817;
x(3)=-0.93827455200273275852364900170872;
x(4)=-0.88641552700440103421315434198220;
x(5)=-0.82000198597390292195394987266975;
x(6)=-0.74012419157855436424382810309998;
x(7)=-0.64809365193697556925249578691075;
x(8)=-0.54542147138883953565837561721837;
x(9)=-0.43379350762604513848708423191335;
x(10)=-0.31504267969616337438679329131981;
x(11)=-0.19111886747361630915863982075707;
x(12)=-0.064056892862605626085043082624745;
x(13)=0.064056892862605626085043082624745;
x(14)=0.19111886747361630915863982075707;
x(15)=0.31504267969616337438679329131981;
x(16)=0.43379350762604513848708423191335;
x(17)=0.54542147138883953565837561721837;
x(18)=0.64809365193697556925249578691075;
x(19)=0.74012419157855436424382810309998;
x(20)=0.82000198597390292195394987266975;
x(21)=0.88641552700440103421315434198220;
x(22)=0.93827455200273275852364900170872;
x(23)=0.97472855597130949819839199300817;
x(24)=0.99518721999702136017999740970074;
elseif (n == 25)
x(1)=-0.99555696979049809790878494689390;
x(2)=-0.97666392145951751149831538647959;
x(3)=-0.94297457122897433941401116965847;
x(4)=-0.89499199787827536885104200678280;
x(5)=-0.83344262876083400142102110869357;
x(6)=-0.75925926303735763057728286520436;
x(7)=-0.67356636847346836448512063324762;
x(8)=-0.57766293024122296772368984161265;
x(9)=-0.47300273144571496052218211500919;
x(10)=-0.36117230580938783773582173012764;
x(11)=-0.24386688372098843204519036279745;
x(12)=-0.12286469261071039638735981880804;
x(13)=0;
x(14)=0.12286469261071039638735981880804;
x(15)=0.24386688372098843204519036279745;
x(16)=0.36117230580938783773582173012764;
x(17)=0.47300273144571496052218211500919;
x(18)=0.57766293024122296772368984161265;
x(19)=0.67356636847346836448512063324762;
x(20)=0.75925926303735763057728286520436;
x(21)=0.83344262876083400142102110869357;
x(22)=0.89499199787827536885104200678280;
x(23)=0.94297457122897433941401116965847;
x(24)=0.97666392145951751149831538647959;
x(25)=0.99555696979049809790878494689390;
elseif (n == 26)
x(1)=-0.99588570114561692900321695932291;
x(2)=-0.97838544595647099110058035431193;
x(3)=-0.94715906666171425013591528351804;
x(4)=-0.90263786198430707421766559923121;
x(5)=-0.84544594278849801879750706146784;
x(6)=-0.77638594882067885619296724724228;
x(7)=-0.69642726041995726486381391372942;
x(8)=-0.60669229301761806323197874691689;
x(9)=-0.50844071482450571769570306472557;
x(10)=-0.40305175512348630648107737709888;
x(11)=-0.29200483948595689514283538207783;
x(12)=-0.17685882035689018396905774841834;
x(13)=-0.059230093429313207093718575198403;
x(14)=0.059230093429313207093718575198403;
x(15)=0.17685882035689018396905774841834;
x(16)=0.29200483948595689514283538207783;
x(17)=0.40305175512348630648107737709888;
x(18)=0.50844071482450571769570306472557;
x(19)=0.60669229301761806323197874691689;
x(20)=0.69642726041995726486381391372942;
x(21)=0.77638594882067885619296724724228;
x(22)=0.84544594278849801879750706146784;
x(23)=0.90263786198430707421766559923121;
x(24)=0.94715906666171425013591528351804;
x(25)=0.97838544595647099110058035431193;
x(26)=0.99588570114561692900321695932291;
elseif (n == 27)
x(1)=-0.99617926288898856693888720838123;
x(2)=-0.97992347596150122285587335566105;
x(3)=-0.95090055781470500685190803064388;
x(4)=-0.90948232067749110430064501820969;
x(5)=-0.85620790801829449030273722270684;
x(6)=-0.79177163907050822714439734410727;
x(7)=-0.71701347373942369929481621164234;
x(8)=-0.63290797194649514092773463763449;
x(9)=-0.54055156457945689490030094155998;
x(10)=-0.44114825175002688058597415568928;
x(11)=-0.33599390363850889973031903420728;
x(12)=-0.22645936543953685885723910736023;
x(13)=-0.11397258560952996693289498386998;
x(14)=0;
x(15)=0.11397258560952996693289498386998;
x(16)=0.22645936543953685885723910736023;
x(17)=0.33599390363850889973031903420728;
x(18)=0.44114825175002688058597415568928;
x(19)=0.54055156457945689490030094155998;
x(20)=0.63290797194649514092773463763449;
x(21)=0.71701347373942369929481621164234;
x(22)=0.79177163907050822714439734410727;
x(23)=0.85620790801829449030273722270684;
x(24)=0.90948232067749110430064501820969;
x(25)=0.95090055781470500685190803064388;
x(26)=0.97992347596150122285587335566105;
x(27)=0.99617926288898856693888720838123;
elseif (n == 28)
x(1)=-0.99644249757395444995043639048331;
x(2)=-0.98130316537087275369455994580783;
x(3)=-0.95425928062893819725410183970522;
x(4)=-0.91563302639213207386968942332993;
x(5)=-0.86589252257439504894225456737969;
x(6)=-0.80564137091717917144788595542528;
x(7)=-0.73561087801363177202814451029253;
x(8)=-0.65665109403886496121989817650674;
x(9)=-0.56972047181140171930800328335643;
x(10)=-0.47587422495511826103441184766743;
x(11)=-0.37625151608907871022135720955609;
x(12)=-0.27206162763517807767682635612577;
x(13)=-0.16456928213338077128147177789117;
x(14)=-0.055079289884034270426516527341880;
x(15)=0.055079289884034270426516527341880;
x(16)=0.16456928213338077128147177789117;
x(17)=0.27206162763517807767682635612577;
x(18)=0.37625151608907871022135720955609;
x(19)=0.47587422495511826103441184766743;
x(20)=0.56972047181140171930800328335643;
x(21)=0.65665109403886496121989817650674;
x(22)=0.73561087801363177202814451029253;
x(23)=0.80564137091717917144788595542528;
x(24)=0.86589252257439504894225456737969;
x(25)=0.91563302639213207386968942332993;
x(26)=0.95425928062893819725410183970522;
x(27)=0.98130316537087275369455994580783;
x(28)=0.99644249757395444995043639048331;
elseif (n == 29)
x(1)=-0.99667944226059658616319153254935;
x(2)=-0.98254550526141317487092601578638;
x(3)=-0.95728559577808772579820803698082;
x(4)=-0.92118023295305878509375343608311;
x(5)=-0.87463780492010279041779342125658;
x(6)=-0.81818548761525244498957221457878;
x(7)=-0.75246285173447713391261007721214;
x(8)=-0.67821453760268651515618500539199;
x(9)=-0.59628179713822782037958621118899;
x(10)=-0.50759295512422764210262791962752;
x(11)=-0.41315288817400866389070658603162;
x(12)=-0.31403163786763993494819592319105;
x(13)=-0.21135228616600107450637572890294;
x(14)=-0.10627823013267923017098239243038;
x(15)=0;
x(16)=0.10627823013267923017098239243038;
x(17)=0.21135228616600107450637572890294;
x(18)=0.31403163786763993494819592319105;
x(19)=0.41315288817400866389070658603162;
x(20)=0.50759295512422764210262791962752;
x(21)=0.59628179713822782037958621118899;
x(22)=0.67821453760268651515618500539199;
x(23)=0.75246285173447713391261007721214;
x(24)=0.81818548761525244498957221457878;
x(25)=0.87463780492010279041779342125658;
x(26)=0.92118023295305878509375343608311;
x(27)=0.95728559577808772579820803698082;
x(28)=0.98254550526141317487092601578638;
x(29)=0.99667944226059658616319153254935;
elseif (n == 30)
x(1)=-0.99689348407464954027163005091870;
x(2)=-0.98366812327974720997003258160566;
x(3)=-0.96002186496830751221687102558180;
x(4)=-0.92620004742927432587932427708047;
x(5)=-0.88256053579205268154311646253023;
x(6)=-0.82956576238276839744289811973250;
x(7)=-0.76777743210482619491797734097450;
x(8)=-0.69785049479331579693229238802664;
x(9)=-0.62052618298924286114047755643119;
x(10)=-0.53662414814201989926416979331107;
x(11)=-0.44703376953808917678060990032285;
x(12)=-0.35270472553087811347103720708937;
x(13)=-0.25463692616788984643980512981781;
x(14)=-0.15386991360858354696379467274326;
x(15)=-0.051471842555317695833025213166723;
x(16)=0.051471842555317695833025213166723;
x(17)=0.15386991360858354696379467274326;
x(18)=0.25463692616788984643980512981781;
x(19)=0.35270472553087811347103720708937;
x(20)=0.44703376953808917678060990032285;
x(21)=0.53662414814201989926416979331107;
x(22)=0.62052618298924286114047755643119;
x(23)=0.69785049479331579693229238802664;
x(24)=0.76777743210482619491797734097450;
x(25)=0.82956576238276839744289811973250;
x(26)=0.88256053579205268154311646253023;
x(27)=0.92620004742927432587932427708047;
x(28)=0.96002186496830751221687102558180;
x(29)=0.98366812327974720997003258160566;
x(30)=0.99689348407464954027163005091870;
elseif (n == 48)
x(1) = -0.9987710072524256;
x(2) = -0.9935301722663504;
x(3) = -0.9841245837228260;
x(4) = -0.9705915925462471;
x(5) = -0.9529877031604301;
x(6) = -0.9313866907065540;
x(7) = -0.9058791367155696;
x(8) = -0.8765720202742477;
x(9) = -0.8435882616243932;
x(10) = -0.8070662040294426;
x(11) = -0.7671590325157405;
x(12) = -0.7240341309238147;
x(13) = -0.6778723796326640;
x(14) = -0.6288673967765133;
x(15) = -0.5772247260839724;
x(16) = -0.5231609747222329;
x(17) = -0.4669029047509581;
x(18) = -0.4086864819907169;
x(19) = -0.3487558862921609;
x(20) = -0.2873624873554553;
x(21) = -0.2247637903946891;
x(22) = -0.1612223560688917;
x(23) = -0.0970046992094628;
x(24) = -0.0323801709628695;
x(25) = 0.0323801709628694;
x(26) = 0.0970046992094629;
x(27) = 0.1612223560688917;
x(28) = 0.2247637903946892;
x(29) = 0.2873624873554556;
x(30) = 0.3487558862921611;
x(31) = 0.4086864819907164;
x(32) = 0.4669029047509585;
x(33) = 0.5231609747222328;
x(34) = 0.5772247260839727;
x(35) = 0.6288673967765133;
x(36) = 0.6778723796326639;
x(37) = 0.7240341309238144;
x(38) = 0.7671590325157402;
x(39) = 0.8070662040294421;
x(40) = 0.8435882616243933;
x(41) = 0.8765720202742481;
x(42) = 0.9058791367155694;
x(43) = 0.9313866907065542;
x(44) = 0.9529877031604304;
x(45) = 0.9705915925462472;
x(46) = 0.9841245837228265;
x(47) = 0.9935301722663505;
x(48) = 0.9987710072524260;
else
error('GAUSS_WEIGHTS - Fatal error! Illegal value of n.');
end
return
end
function w = gauss_weights(n);
%
% function w = gauss_weights(n);
%
% For 1 <= n <= 30, returns the weights w of an
% n point Gauss-Legendre quadrature rule over the interval [-1,1].
%
w = ones(1,n);
if ( n == 1 )
w(1) = 2.0;
elseif ( n == 2 )
w(1) = 1.0;
w(2) = w(1);
elseif ( n == 3 )
w(1) = 0.5555555555555555555555555555565;
w(2) = 0.8888888888888888888888888888889;
w(3) = 0.5555555555555555555555555555565;
elseif ( n == 4 )
w(1) = 0.347854845137453857373063949222;
w(2) = 0.652145154862546142626936050778;
w(3) = 0.652145154862546142626936050778;
w(4) = 0.347854845137453857373063949222;
elseif ( n == 5 )
w(1) = 0.236926885056189087514264040720;
w(2) = 0.478628670499366468041291514836;
w(3) = 0.568888888888888888888888888889;
w(4) = 0.478628670499366468041291514836;
w(5) = 0.236926885056189087514264040720;
elseif ( n == 6 )
w(1) = 0.171324492379170345040296142173;
w(2) = 0.360761573048138607569833513838;
w(3) = 0.467913934572691047389870343990;
w(4) = 0.467913934572691047389870343990;
w(5) = 0.360761573048138607569833513838;
w(6) = 0.171324492379170345040296142173;
elseif ( n == 7 )
w(1) = 0.129484966168869693270611432679;
w(2) = 0.279705391489276667901467771424;
w(3) = 0.381830050505118944950369775489;
w(4) = 0.417959183673469387755102040816;
w(5) = 0.381830050505118944950369775489;
w(6) = 0.279705391489276667901467771424;
w(7) = 0.129484966168869693270611432679;
elseif ( n == 8 )
w(1) = 0.101228536290376259152531354310;
w(2) = 0.222381034453374470544355994426;
w(3) = 0.313706645877887287337962201987;
w(4) = 0.362683783378361982965150449277;
w(5) = 0.362683783378361982965150449277;
w(6) = 0.313706645877887287337962201987;
w(7) = 0.222381034453374470544355994426;
w(8) = 0.101228536290376259152531354310;
elseif ( n == 9 )
w(1) = 0.0812743883615744119718921581105;
w(2) = 0.180648160694857404058472031243;
w(3) = 0.260610696402935462318742869419;
w(4) = 0.312347077040002840068630406584;
w(5) = 0.330239355001259763164525069287;
w(6) = 0.312347077040002840068630406584;
w(7) = 0.260610696402935462318742869419;
w(8) = 0.180648160694857404058472031243;
w(9) = 0.0812743883615744119718921581105;
elseif ( n == 10 )
w(1) = 0.0666713443086881375935688098933;
w(2) = 0.149451349150580593145776339658;
w(3) = 0.219086362515982043995534934228;
w(4) = 0.269266719309996355091226921569;
w(5) = 0.295524224714752870173892994651;
w(6) = 0.295524224714752870173892994651;
w(7) = 0.269266719309996355091226921569;
w(8) = 0.219086362515982043995534934228;
w(9) = 0.149451349150580593145776339658;
w(10) = 0.0666713443086881375935688098933;
elseif ( n == 11 )
w(1)= 0.0556685671161736664827537204425;
w(2)= 0.125580369464904624634694299224;
w(3)= 0.186290210927734251426097641432;
w(4)= 0.233193764591990479918523704843;
w(5)= 0.262804544510246662180688869891;
w(6)= 0.272925086777900630714483528336;
w(7)= 0.262804544510246662180688869891;
w(8)= 0.233193764591990479918523704843;
w(9)= 0.186290210927734251426097641432;
w(10)= 0.125580369464904624634694299224;
w(11)= 0.0556685671161736664827537204425;
elseif ( n == 12 )
w(1)= 0.0471753363865118271946159614850;
w(2)= 0.106939325995318430960254718194;
w(3)= 0.160078328543346226334652529543;
w(4)= 0.203167426723065921749064455810;
w(5)= 0.233492536538354808760849898925;
w(6)= 0.249147045813402785000562436043;
w(7)= 0.249147045813402785000562436043;
w(8)= 0.233492536538354808760849898925;
w(9)= 0.203167426723065921749064455810;
w(10)= 0.160078328543346226334652529543;
w(11)= 0.106939325995318430960254718194;
w(12)= 0.0471753363865118271946159614850;
elseif ( n == 13 )
w(1)= 0.0404840047653158795200215922010;
w(2)= 0.0921214998377284479144217759538;
w(3)= 0.138873510219787238463601776869;
w(4)= 0.178145980761945738280046691996;
w(5)= 0.207816047536888502312523219306;
w(6)= 0.226283180262897238412090186040;
w(7)= 0.232551553230873910194589515269;
w(8)= 0.226283180262897238412090186040;
w(9)= 0.207816047536888502312523219306;
w(10)= 0.178145980761945738280046691996;
w(11)= 0.138873510219787238463601776869;
w(12)= 0.0921214998377284479144217759538;
w(13)= 0.0404840047653158795200215922010;
elseif ( n == 14 )
w(1)= 0.0351194603317518630318328761382;
w(2)= 0.0801580871597602098056332770629;
w(3)= 0.121518570687903184689414809072;
w(4)= 0.157203167158193534569601938624;
w(5)= 0.185538397477937813741716590125;
w(6)= 0.205198463721295603965924065661;
w(7)= 0.215263853463157790195876443316;
w(8)= 0.215263853463157790195876443316;
w(9)= 0.205198463721295603965924065661;
w(10)= 0.185538397477937813741716590125;
w(11)= 0.157203167158193534569601938624;
w(12)= 0.121518570687903184689414809072;
w(13)= 0.0801580871597602098056332770629;
w(14)= 0.0351194603317518630318328761382;
elseif ( n == 15 )
w(1)= 0.0307532419961172683546283935772;
w(2)= 0.0703660474881081247092674164507;
w(3)= 0.107159220467171935011869546686;
w(4)= 0.139570677926154314447804794511;
w(5)= 0.166269205816993933553200860481;
w(6)= 0.186161000015562211026800561866;
w(7)= 0.198431485327111576456118326444;
w(8)= 0.202578241925561272880620199968;
w(9)= 0.198431485327111576456118326444;
w(10)= 0.186161000015562211026800561866;
w(11)= 0.166269205816993933553200860481;
w(12)= 0.139570677926154314447804794511;
w(13)= 0.107159220467171935011869546686;
w(14)= 0.0703660474881081247092674164507;
w(15)= 0.0307532419961172683546283935772;
elseif ( n == 16 )
w(1)= 0.0271524594117540948517805724560;
w(2)= 0.0622535239386478928628438369944;
w(3)= 0.0951585116824927848099251076022;
w(4)= 0.124628971255533872052476282192;
w(5)= 0.149595988816576732081501730547;
w(6)= 0.169156519395002538189312079030;
w(7)= 0.182603415044923588866763667969;
w(8)= 0.189450610455068496285396723208;
w(9)= 0.189450610455068496285396723208;
w(10)= 0.182603415044923588866763667969;
w(11)= 0.169156519395002538189312079030;
w(12)= 0.149595988816576732081501730547;
w(13)= 0.124628971255533872052476282192;
w(14)= 0.0951585116824927848099251076022;
w(15)= 0.0622535239386478928628438369944;
w(16)= 0.0271524594117540948517805724560;
elseif ( n == 17 )
w(1)= 0.0241483028685479319601100262876;
w(2)= 0.0554595293739872011294401653582;
w(3)= 0.0850361483171791808835353701911;
w(4)= 0.111883847193403971094788385626;
w(5)= 0.135136368468525473286319981702;
w(6)= 0.154045761076810288081431594802;
w(7)= 0.168004102156450044509970663788;
w(8)= 0.176562705366992646325270990113;
w(9)= 0.179446470356206525458265644262;
w(10)= 0.176562705366992646325270990113;
w(11)= 0.168004102156450044509970663788;
w(12)= 0.154045761076810288081431594802;
w(13)= 0.135136368468525473286319981702;
w(14)= 0.111883847193403971094788385626;
w(15)= 0.0850361483171791808835353701911;
w(16)= 0.0554595293739872011294401653582;
w(17)= 0.0241483028685479319601100262876;
elseif ( n == 18 )
w(1)= 0.0216160135264833103133427102665;
w(2)= 0.0497145488949697964533349462026;
w(3)= 0.0764257302548890565291296776166;
w(4)= 0.100942044106287165562813984925;
w(5)= 0.122555206711478460184519126800;
w(6)= 0.140642914670650651204731303752;
w(7)= 0.154684675126265244925418003836;
w(8)= 0.164276483745832722986053776466;
w(9)= 0.169142382963143591840656470135;
w(10)= 0.169142382963143591840656470135;
w(11)= 0.164276483745832722986053776466;
w(12)= 0.154684675126265244925418003836;
w(13)= 0.140642914670650651204731303752;
w(14)= 0.122555206711478460184519126800;
w(15)= 0.100942044106287165562813984925;
w(16)= 0.0764257302548890565291296776166;
w(17)= 0.0497145488949697964533349462026;
w(18)= 0.0216160135264833103133427102665;
elseif ( n == 19 )
w(1)= 0.0194617882297264770363120414644;
w(2)= 0.0448142267656996003328381574020;
w(3)= 0.0690445427376412265807082580060;
w(4)= 0.0914900216224499994644620941238;
w(5)= 0.111566645547333994716023901682;
w(6)= 0.128753962539336227675515784857;
w(7)= 0.142606702173606611775746109442;
w(8)= 0.152766042065859666778855400898;
w(9)= 0.158968843393954347649956439465;
w(10)= 0.161054449848783695979163625321;
w(11)= 0.158968843393954347649956439465;
w(12)= 0.152766042065859666778855400898;
w(13)= 0.142606702173606611775746109442;
w(14)= 0.128753962539336227675515784857;
w(15)= 0.111566645547333994716023901682;
w(16)= 0.0914900216224499994644620941238;
w(17)= 0.0690445427376412265807082580060;
w(18)= 0.0448142267656996003328381574020;
w(19)= 0.0194617882297264770363120414644;
elseif ( n == 20 )
w(1)= 0.0176140071391521183118619623519;
w(2)= 0.0406014298003869413310399522749;
w(3)= 0.0626720483341090635695065351870;
w(4)= 0.0832767415767047487247581432220;
w(5)= 0.101930119817240435036750135480;
w(6)= 0.118194531961518417312377377711;
w(7)= 0.131688638449176626898494499748;
w(8)= 0.142096109318382051329298325067;
w(9)= 0.149172986472603746787828737002;
w(10)= 0.152753387130725850698084331955;
w(11)= 0.152753387130725850698084331955;
w(12)= 0.149172986472603746787828737002;
w(13)= 0.142096109318382051329298325067;
w(14)= 0.131688638449176626898494499748;
w(15)= 0.118194531961518417312377377711;
w(16)= 0.101930119817240435036750135480;
w(17)= 0.0832767415767047487247581432220;
w(18)= 0.0626720483341090635695065351870;
w(19)= 0.0406014298003869413310399522749;
w(20)= 0.0176140071391521183118619623519;
elseif ( n == 21 )
w(1)=0.016017228257774333324224616858471;
w(2)=0.036953789770852493799950668299330;
w(3)=0.057134425426857208283635826472448;
w(4)=0.076100113628379302017051653300183;
w(5)=0.093444423456033861553289741113932;
w(6)=0.10879729916714837766347457807011;
w(7)=0.12183141605372853419536717712573;
w(8)=0.13226893863333746178105257449678;
w(9)=0.13988739479107315472213342386758;
w(10)=0.14452440398997005906382716655375;
w(11)=0.14608113364969042719198514768337;
w(12)=0.14452440398997005906382716655375;
w(13)=0.13988739479107315472213342386758;
w(14)=0.13226893863333746178105257449678;
w(15)=0.12183141605372853419536717712573;
w(16)=0.10879729916714837766347457807011;
w(17)=0.093444423456033861553289741113932;
w(18)=0.076100113628379302017051653300183;
w(19)=0.057134425426857208283635826472448;
w(20)=0.036953789770852493799950668299330;
w(21)=0.016017228257774333324224616858471;
elseif ( n == 22 )
w(1)=0.014627995298272200684991098047185;
w(2)=0.033774901584814154793302246865913;
w(3)=0.052293335152683285940312051273211;
w(4)=0.069796468424520488094961418930218;
w(5)=0.085941606217067727414443681372703;
w(6)=0.10041414444288096493207883783054;
w(7)=0.11293229608053921839340060742178;
w(8)=0.12325237681051242428556098615481;
w(9)=0.13117350478706237073296499253031;
w(10)=0.13654149834601517135257383123152;
w(11)=0.13925187285563199337541024834181;
w(12)=0.13925187285563199337541024834181;
w(13)=0.13654149834601517135257383123152;
w(14)=0.13117350478706237073296499253031;
w(15)=0.12325237681051242428556098615481;
w(16)=0.11293229608053921839340060742178;
w(17)=0.10041414444288096493207883783054;
w(18)=0.085941606217067727414443681372703;
w(19)=0.069796468424520488094961418930218;
w(20)=0.052293335152683285940312051273211;
w(21)=0.033774901584814154793302246865913;
w(22)=0.014627995298272200684991098047185;
elseif ( n == 23 )
w(1)=0.013411859487141772081309493458615;
w(2)=0.030988005856979444310694219641885;
w(3)=0.048037671731084668571641071632034;
w(4)=0.064232421408525852127169615158911;
w(5)=0.079281411776718954922892524742043;
w(6)=0.092915766060035147477018617369765;
w(7)=0.10489209146454141007408618501474;
w(8)=0.11499664022241136494164351293396;
w(9)=0.12304908430672953046757840067201;
w(10)=0.12890572218808214997859533939979;
w(11)=0.13246203940469661737164246470332;
w(12)=0.13365457218610617535145711054584;
w(13)=0.13246203940469661737164246470332;
w(14)=0.12890572218808214997859533939979;
w(15)=0.12304908430672953046757840067201;
w(16)=0.11499664022241136494164351293396;
w(17)=0.10489209146454141007408618501474;
w(18)=0.092915766060035147477018617369765;
w(19)=0.079281411776718954922892524742043;
w(20)=0.064232421408525852127169615158911;
w(21)=0.048037671731084668571641071632034;
w(22)=0.030988005856979444310694219641885;
w(23)=0.013411859487141772081309493458615;
elseif ( n == 24 )
w(1)=0.012341229799987199546805667070037;
w(2)=0.028531388628933663181307815951878;
w(3)=0.044277438817419806168602748211338;
w(4)=0.059298584915436780746367758500109;
w(5)=0.073346481411080305734033615253117;
w(6)=0.086190161531953275917185202983743;
w(7)=0.097618652104113888269880664464247;
w(8)=0.10744427011596563478257734244661;
w(9)=0.11550566805372560135334448390678;
w(10)=0.12167047292780339120446315347626;
w(11)=0.12583745634682829612137538251118;
w(12)=0.12793819534675215697405616522470;
w(13)=0.12793819534675215697405616522470;
w(14)=0.12583745634682829612137538251118;
w(15)=0.12167047292780339120446315347626;
w(16)=0.11550566805372560135334448390678;
w(17)=0.10744427011596563478257734244661;
w(18)=0.097618652104113888269880664464247;
w(19)=0.086190161531953275917185202983743;
w(20)=0.073346481411080305734033615253117;
w(21)=0.059298584915436780746367758500109;
w(22)=0.044277438817419806168602748211338;
w(23)=0.028531388628933663181307815951878;
w(24)=0.012341229799987199546805667070037;
elseif ( n == 25 )
w(1)=0.011393798501026287947902964113235;
w(2)=0.026354986615032137261901815295299;
w(3)=0.040939156701306312655623487711646;
w(4)=0.054904695975835191925936891540473;
w(5)=0.068038333812356917207187185656708;
w(6)=0.080140700335001018013234959669111;
w(7)=0.091028261982963649811497220702892;
w(8)=0.10053594906705064420220689039269;
w(9)=0.10851962447426365311609395705012;
w(10)=0.11485825914571164833932554586956;
w(11)=0.11945576353578477222817812651290;
w(12)=0.12224244299031004168895951894585;
w(13)=0.12317605372671545120390287307905;
w(14)=0.12224244299031004168895951894585;
w(15)=0.11945576353578477222817812651290;
w(16)=0.11485825914571164833932554586956;
w(17)=0.10851962447426365311609395705012;
w(18)=0.10053594906705064420220689039269;
w(19)=0.091028261982963649811497220702892;
w(20)=0.080140700335001018013234959669111;
w(21)=0.068038333812356917207187185656708;
w(22)=0.054904695975835191925936891540473;
w(23)=0.040939156701306312655623487711646;
w(24)=0.026354986615032137261901815295299;
w(25)=0.011393798501026287947902964113235;
elseif ( n == 26 )
w(1)=0.010551372617343007155651187685252;
w(2)=0.024417851092631908789615827519788;
w(3)=0.037962383294362763950303141248850;
w(4)=0.050975825297147811998319900724073;
w(5)=0.063274046329574835539453689907045;
w(6)=0.074684149765659745887075796102848;
w(7)=0.085045894313485239210447765079982;
w(8)=0.094213800355914148463664883067303;
w(9)=0.10205916109442542323841407025343;
w(10)=0.10847184052857659065657942672800;
w(11)=0.11336181654631966654944071844260;
w(12)=0.11666044348529658204466250754036;
w(13)=0.11832141527926227651637108570047;
w(14)=0.11832141527926227651637108570047;
w(15)=0.11666044348529658204466250754036;
w(16)=0.11336181654631966654944071844260;
w(17)=0.10847184052857659065657942672800;
w(18)=0.10205916109442542323841407025343;
w(19)=0.094213800355914148463664883067303;
w(20)=0.085045894313485239210447765079982;
w(21)=0.074684149765659745887075796102848;
w(22)=0.063274046329574835539453689907045;
w(23)=0.050975825297147811998319900724073;
w(24)=0.037962383294362763950303141248850;
w(25)=0.024417851092631908789615827519788;
w(26)=0.010551372617343007155651187685252;
elseif ( n == 27 )
w(1)=0.0097989960512943602611500550912591;
w(2)=0.022686231596180623196034206446761;
w(3)=0.035297053757419711022578289304712;
w(4)=0.047449412520615062704096710114185;
w(5)=0.058983536859833599110300833719532;
w(6)=0.069748823766245592984322888356667;
w(7)=0.079604867773057771263074959009842;
w(8)=0.088423158543756950194322802853749;
w(9)=0.096088727370028507565652646558106;
w(10)=0.10250163781774579867124771153266;
w(11)=0.10757828578853318721216298442666;
w(12)=0.11125248835684519267216309604285;
w(13)=0.11347634610896514862036994809210;
w(14)=0.11422086737895698904504573690184;
w(15)=0.11347634610896514862036994809210;
w(16)=0.11125248835684519267216309604285;
w(17)=0.10757828578853318721216298442666;
w(18)=0.10250163781774579867124771153266;
w(19)=0.096088727370028507565652646558106;
w(20)=0.088423158543756950194322802853749;
w(21)=0.079604867773057771263074959009842;
w(22)=0.069748823766245592984322888356667;
w(23)=0.058983536859833599110300833719532;
w(24)=0.047449412520615062704096710114185;
w(25)=0.035297053757419711022578289304712;
w(26)=0.022686231596180623196034206446761;
w(27)=0.0097989960512943602611500550912591;
elseif ( n == 28 )
w(1)=0.0091242825930945177388161539229517;
w(2)=0.021132112592771259751500380993265;
w(3)=0.032901427782304379977630819170532;
w(4)=0.044272934759004227839587877653207;
w(5)=0.055107345675716745431482918226946;
w(6)=0.065272923966999595793397566775505;
w(7)=0.074646214234568779023931887173022;
w(8)=0.083113417228901218390396498244332;
w(9)=0.090571744393032840942186031336784;
w(10)=0.096930657997929915850489006095441;
w(11)=0.10211296757806076981421663850571;
w(12)=0.10605576592284641791041643699681;
w(13)=0.10871119225829413525357151930367;
w(14)=0.11004701301647519628237626560182;
w(15)=0.11004701301647519628237626560182;
w(16)=0.10871119225829413525357151930367;
w(17)=0.10605576592284641791041643699681;
w(18)=0.10211296757806076981421663850571;
w(19)=0.096930657997929915850489006095441;
w(20)=0.090571744393032840942186031336784;
w(21)=0.083113417228901218390396498244332;
w(22)=0.074646214234568779023931887173022;
w(23)=0.065272923966999595793397566775505;
w(24)=0.055107345675716745431482918226946;
w(25)=0.044272934759004227839587877653207;
w(26)=0.032901427782304379977630819170532;
w(27)=0.021132112592771259751500380993265;
w(28)=0.0091242825930945177388161539229517;
elseif ( n == 29 )
w(1)=0.0085169038787464096542638133022498;
w(2)=0.019732085056122705983859801640396;
w(3)=0.030740492202093622644408525374617;
w(4)=0.041402062518682836104830010114077;
w(5)=0.051594826902497923912594381179543;
w(6)=0.061203090657079138542109848023907;
w(7)=0.070117933255051278569581486948879;
w(8)=0.078238327135763783828144888659680;
w(9)=0.085472257366172527545344849297208;
w(10)=0.091737757139258763347966411077111;
w(11)=0.096963834094408606301900074882689;
w(12)=0.10109127375991496612182054690750;
w(13)=0.10407331007772937391332847128512;
w(14)=0.10587615509732094140659132785219;
w(15)=0.10647938171831424424651112690968;
w(16)=0.10587615509732094140659132785219;
w(17)=0.10407331007772937391332847128512;
w(18)=0.10109127375991496612182054690750;
w(19)=0.096963834094408606301900074882689;
w(20)=0.091737757139258763347966411077111;
w(21)=0.085472257366172527545344849297208;
w(22)=0.078238327135763783828144888659680;
w(23)=0.070117933255051278569581486948879;
w(24)=0.061203090657079138542109848023907;
w(25)=0.051594826902497923912594381179543;
w(26)=0.041402062518682836104830010114077;
w(27)=0.030740492202093622644408525374617;
w(28)=0.019732085056122705983859801640396;
w(29)=0.0085169038787464096542638133022498;
elseif ( n == 30 )
w(1)=0.0079681924961666056154658834746736;
w(2)=0.018466468311090959142302131912047;
w(3)=0.028784707883323369349719179611292;
w(4)=0.038799192569627049596801936446348;
w(5)=0.048402672830594052902938140422808;
w(6)=0.057493156217619066481721689402056;
w(7)=0.065974229882180495128128515115962;
w(8)=0.073755974737705206268243850022191;
w(9)=0.080755895229420215354694938460530;
w(10)=0.086899787201082979802387530715126;
w(11)=0.092122522237786128717632707087619;
w(12)=0.096368737174644259639468626351810;
w(13)=0.099593420586795267062780282103569;
w(14)=0.10176238974840550459642895216855;
w(15)=0.10285265289355884034128563670542;
w(16)=0.10285265289355884034128563670542;
w(17)=0.10176238974840550459642895216855;
w(18)=0.099593420586795267062780282103569;
w(19)=0.096368737174644259639468626351810;
w(20)=0.092122522237786128717632707087619;
w(21)=0.086899787201082979802387530715126;
w(22)=0.080755895229420215354694938460530;
w(23)=0.073755974737705206268243850022191;
w(24)=0.065974229882180495128128515115962;
w(25)=0.057493156217619066481721689402056;
w(26)=0.048402672830594052902938140422808;
w(27)=0.038799192569627049596801936446348;
w(28)=0.028784707883323369349719179611292;
w(29)=0.018466468311090959142302131912047;
w(30)=0.0079681924961666056154658834746736;
elseif (n == 48)
w(1) = 0.0031533460523059;
w(2) = 0.0073275539012765;
w(3) = 0.0114772345792347;
w(4) = 0.0155793157229437;
w(5) = 0.0196161604573558;
w(6) = 0.0235707608393237;
w(7) = 0.0274265097083563;
w(8) = 0.0311672278327980;
w(9) = 0.0347772225647707;
w(10) = 0.0382413510658304;
w(11) = 0.0415450829434656;
w(12) = 0.0446745608566937;
w(13) = 0.0476166584924910;
w(14) = 0.0503590355538545;
w(15) = 0.0528901894851937;
w(16) = 0.0551995036999840;
w(17) = 0.0572772921004029;
w(18) = 0.0591148396983953;
w(19) = 0.0607044391658945;
w(20) = 0.0620394231598924;
w(21) = 0.0631141922862541;
w(22) = 0.0639242385846489;
w(23) = 0.0644661644359506;
w(24) = 0.0647376968126839;
w(25) = 0.0647376968126836;
w(26) = 0.0644661644359505;
w(27) = 0.0639242385846487;
w(28) = 0.0631141922862533;
w(29) = 0.0620394231598915;
w(30) = 0.0607044391658939;
w(31) = 0.0591148396983964;
w(32) = 0.0572772921004031;
w(33) = 0.0551995036999844;
w(34) = 0.0528901894851943;
w(35) = 0.0503590355538540;
w(36) = 0.0476166584924911;
w(37) = 0.0446745608566940;
w(38) = 0.0415450829434649;
w(39) = 0.0382413510658305;
w(40) = 0.0347772225647698;
w(41) = 0.0311672278327982;
w(42) = 0.0274265097083569;
w(43) = 0.0235707608393248;
w(44) = 0.0196161604573554;
w(45) = 0.0155793157229439;
w(46) = 0.0114772345792341;
w(47) = 0.0073275539012765;
w(48) = 0.0031533460523052;
else
error('GAUSS_WEIGHTS - Fatal error! Illegal value of n.');
end
return
end
function xi = lobatto_points(n)
xi = ones(1,n);
switch n
case(2)
xi(1)=-1.0000000000000000000000000000000;
xi(2)=1.0000000000000000000000000000000;
case(3)
xi(1)=-1.0000000000000000000000000000000;
xi(2)=0.0;
xi(3)=1.0000000000000000000000000000000;
case(4)
xi(1)=-1.0000000000000000000000000000000;
xi(2)=-0.44721359549995793928183473374626;
xi(3)=0.44721359549995793928183473374626;
xi(4)=1.0000000000000000000000000000000;
case(5)
xi(1)=-1.0000000000000000000000000000000;
xi(2)=-0.65465367070797714379829245624686;
xi(3)=0.0;
xi(4)=0.65465367070797714379829245624686;
xi(5)=1.0000000000000000000000000000000;
case(6)
xi(1)=-1.0000000000000000000000000000000;
xi(2)=-0.76505532392946469285100297395934;
xi(3)=-0.28523151648064509631415099404088;
xi(4)=0.28523151648064509631415099404088;
xi(5)=0.76505532392946469285100297395934;
xi(6)=1.0000000000000000000000000000000;
case(7)
xi(1)=-1.0000000000000000000000000000000;
xi(2)=-0.83022389627856692987203221396747;
xi(3)=-0.46884879347071421380377188190877;
xi(4)=0.0;
xi(5)=0.46884879347071421380377188190877;
xi(6)=0.83022389627856692987203221396747;
xi(7)=1.0000000000000000000000000000000;
case(8)
xi(1)=-1.0000000000000000000000000000000;
xi(2)=-0.87174014850960661533744576122066;
xi(3)=-0.59170018143314230214451073139795;
xi(4)=-0.20929921790247886876865726034535;
xi(5)=0.20929921790247886876865726034535;
xi(6)=0.59170018143314230214451073139795;
xi(7)=0.87174014850960661533744576122066;
xi(8)=1.0000000000000000000000000000000;
case(9)
xi(1)=-1.0000000000000000000000000000000;
xi(2)=-0.89975799541146015731234524441834;
xi(3)=-0.67718627951073775344588542709134;
xi(4)=-0.36311746382617815871075206870866;
xi(5)=0.0;
xi(6)=0.36311746382617815871075206870866;
xi(7)=0.67718627951073775344588542709134;
xi(8)=0.89975799541146015731234524441834;
xi(9)=1.0000000000000000000000000000000;
case(10)
xi(1)=-1.0000000000000000000000000000000;
xi(2)=-0.91953390816645881382893266082234;
xi(3)=-0.73877386510550507500310617485983;
xi(4)=-0.47792494981044449566117509273126;
xi(5)=-0.16527895766638702462621976595817;
xi(6)=0.16527895766638702462621976595817;
xi(7)=0.47792494981044449566117509273126;
xi(8)=0.73877386510550507500310617485983;
xi(9)=0.91953390816645881382893266082234;
xi(10)=1.0000000000000000000000000000000;
case(11)
xi(1)=-1.0000000000000000000000000000000;
xi(2)=-0.93400143040805913433227413609938;
xi(3)=-0.78448347366314441862241781610846;
xi(4)=-0.56523532699620500647096396947775;
xi(5)=-0.29575813558693939143191151555906;
xi(6)=0.0;
xi(7)=0.29575813558693939143191151555906;
xi(8)=0.56523532699620500647096396947775;
xi(9)=0.78448347366314441862241781610846;
xi(10)=0.93400143040805913433227413609938;
xi(11)=1.0000000000000000000000000000000;
case(12)
xi(1)=-1.0000000000000000000000000000000;
xi(2)=-0.94489927222288222340758013830322;
xi(3)=-0.81927932164400667834864158171690;
xi(4)=-0.63287615303186067766240485444366;
xi(5)=-0.39953094096534893226434979156697;
xi(6)=-0.13655293285492755486406185573969;
xi(7)=0.13655293285492755486406185573969;
xi(8)=0.39953094096534893226434979156697;
xi(9)=0.63287615303186067766240485444366;
xi(10)=0.81927932164400667834864158171690;
xi(11)=0.94489927222288222340758013830322;
xi(12)=1.0000000000000000000000000000000;
case(13)
xi(1)=-1.0000000000000000000000000000000;
xi(2)=-0.95330984664216391189690546475545;
xi(3)=-0.84634756465187231686592560709875;
xi(4)=-0.68618846908175742607275903956636;
xi(5)=-0.48290982109133620174693723363693;
xi(6)=-0.24928693010623999256867370037423;
xi(7)=0.0;
xi(8)=0.24928693010623999256867370037423;
xi(9)=0.48290982109133620174693723363693;
xi(10)=0.68618846908175742607275903956636;
xi(11)=0.84634756465187231686592560709875;
xi(12)=0.95330984664216391189690546475545;
xi(13)=1.0000000000000000000000000000000;
case(14)
xi(1)=-1.0000000000000000000000000000000;
xi(2)=-0.95993504526726090135510016201542;
xi(3)=-0.86780105383034725100022020290826;
xi(4)=-0.72886859909132614058467240052088;
xi(5)=-0.55063940292864705531662270585908;
xi(6)=-0.34272401334271284504390340364167;
xi(7)=-0.11633186888370386765877670973616;
xi(8)=0.11633186888370386765877670973616;
xi(9)=0.34272401334271284504390340364167;
xi(10)=0.55063940292864705531662270585908;
xi(11)=0.72886859909132614058467240052088;
xi(12)=0.86780105383034725100022020290826;
xi(13)=0.95993504526726090135510016201542;
xi(14)=1.0000000000000000000000000000000;
case(15)
xi(1)=-1.0000000000000000000000000000000;
xi(2)=-0.96524592650383857279585139206960;
xi(3)=-0.88508204422297629882540163148223;
xi(4)=-0.76351968995181520070411847597629;
xi(5)=-0.60625320546984571112352993863673;
xi(6)=-0.42063805471367248092189693873858;
xi(7)=-0.21535395536379423822567944627292;
xi(8)=0.0;
xi(9)=0.21535395536379423822567944627292;
xi(10)=0.42063805471367248092189693873858;
xi(11)=0.60625320546984571112352993863673;
xi(12)=0.76351968995181520070411847597629;
xi(13)=0.88508204422297629882540163148223;
xi(14)=0.96524592650383857279585139206960;
xi(15)=1.0000000000000000000000000000000;
case(16)
xi(1)=-1.0000000000000000000000000000000;
xi(2)=-0.96956804627021793295224273836746;
xi(3)=-0.89920053309347209299462826151985;
xi(4)=-0.79200829186181506393108827096315;
xi(5)=-0.65238870288249308946788321964058;
xi(6)=-0.48605942188713761178189078584687;
xi(7)=-0.29983046890076320809835345472230;
xi(8)=-0.10132627352194944784303300504592;
xi(9)=0.10132627352194944784303300504592;
xi(10)=0.29983046890076320809835345472230;
xi(11)=0.48605942188713761178189078584687;
xi(12)=0.65238870288249308946788321964058;
xi(13)=0.79200829186181506393108827096315;
xi(14)=0.89920053309347209299462826151985;
xi(15)=0.96956804627021793295224273836746;
xi(16)=1.0000000000000000000000000000000;
case(17)
xi(1)=-1.0000000000000000000000000000000;
xi(2)=-0.97313217663141831415697950187372;
xi(3)=-0.91087999591557359562380250639773;
xi(4)=-0.81569625122177030710675055323753;
xi(5)=-0.69102898062768470539491935737245;
xi(6)=-0.54138539933010153912373340750406;
xi(7)=-0.37217443356547704190723468073526;
xi(8)=-0.18951197351831738830426301475311;
xi(9)=0.0;
xi(10)=0.18951197351831738830426301475311;
xi(11)=0.37217443356547704190723468073526;
xi(12)=0.54138539933010153912373340750406;
xi(13)=0.69102898062768470539491935737245;
xi(14)=0.81569625122177030710675055323753;
xi(15)=0.91087999591557359562380250639773;
xi(16)=0.97313217663141831415697950187372;
xi(17)=1.0000000000000000000000000000000;
case(18)
xi(1)=-1.0000000000000000000000000000000;
xi(2)=-0.97610555741219854286451892434170;
xi(3)=-0.92064918534753387383785462543128;
xi(4)=-0.83559353521809021371364636232794;
xi(5)=-0.72367932928324268130621036530207;
xi(6)=-0.58850483431866176117353589319356;
xi(7)=-0.43441503691212397534228713674067;
xi(8)=-0.26636265287828098416766533202560;
xi(9)=-0.089749093484652111022645010088562;
xi(10)=0.089749093484652111022645010088562;
xi(11)=0.26636265287828098416766533202560;
xi(12)=0.43441503691212397534228713674067;
xi(13)=0.58850483431866176117353589319356;
xi(14)=0.72367932928324268130621036530207;
xi(15)=0.83559353521809021371364636232794;
xi(16)=0.92064918534753387383785462543128;
xi(17)=0.97610555741219854286451892434170;
xi(18)=1.0000000000000000000000000000000;
case(19)
xi(1)=-1.0000000000000000000000000000000;
xi(2)=-0.97861176622208009515263406311022;
xi(3)=-0.92890152815258624371794025879655;
xi(4)=-0.85246057779664609308595597004106;
xi(5)=-0.75149420255261301416363748963394;
xi(6)=-0.62890813726522049776683230622873;
xi(7)=-0.48822928568071350277790963762492;
xi(8)=-0.33350484782449861029850010384493;
xi(9)=-0.16918602340928157137515415344488;
xi(10)=0.0;
xi(11)=0.16918602340928157137515415344488;
xi(12)=0.33350484782449861029850010384493;
xi(13)=0.48822928568071350277790963762492;
xi(14)=0.62890813726522049776683230622873;
xi(15)=0.75149420255261301416363748963394;
xi(16)=0.85246057779664609308595597004106;
xi(17)=0.92890152815258624371794025879655;
xi(18)=0.97861176622208009515263406311022;
xi(19)=1.0000000000000000000000000000000;
case(20)
xi(1)=-1.0000000000000000000000000000000;
xi(2)=-0.98074370489391417192544643858423;
xi(3)=-0.93593449881266543571618158493063;
xi(4)=-0.86687797808995014130984721461629;
xi(5)=-0.77536826095205587041431752759469;
xi(6)=-0.66377640229031128984640332297116;
xi(7)=-0.53499286403188626164813596182898;
xi(8)=-0.39235318371390929938647470381582;
xi(9)=-0.23955170592298649518240135692709;
xi(10)=-0.080545937238821837975944518159554;
xi(11)=0.080545937238821837975944518159554;
xi(12)=0.23955170592298649518240135692709;
xi(13)=0.39235318371390929938647470381582;
xi(14)=0.53499286403188626164813596182898;
xi(15)=0.66377640229031128984640332297116;
xi(16)=0.77536826095205587041431752759469;
xi(17)=0.86687797808995014130984721461629;
xi(18)=0.93593449881266543571618158493063;
xi(19)=0.98074370489391417192544643858423;
xi(20)=1.0000000000000000000000000000000;
case(21)
xi(1)=-1.0000000000000000000000000000000;
xi(2)=-0.98257229660454802823448127655541;
xi(3)=-0.94197629695974553429610265066144;
xi(4)=-0.87929475532359046445115359630494;
xi(5)=-0.79600192607771240474431258966036;
xi(6)=-0.69405102606222323262731639319467;
xi(7)=-0.57583196026183068692702187033809;
xi(8)=-0.44411578327900210119451634960735;
xi(9)=-0.30198985650876488727535186785875;
xi(10)=-0.15278551580218546600635832848567;
xi(11)=0.0;
xi(12)=0.15278551580218546600635832848567;
xi(13)=0.30198985650876488727535186785875;
xi(14)=0.44411578327900210119451634960735;
xi(15)=0.57583196026183068692702187033809;
xi(16)=0.69405102606222323262731639319467;
xi(17)=0.79600192607771240474431258966036;
xi(18)=0.87929475532359046445115359630494;
xi(19)=0.94197629695974553429610265066144;
xi(20)=0.98257229660454802823448127655541;
xi(21)=1.0000000000000000000000000000000;
case(22)
xi(1)=-1.0000000000000000000000000000000;
xi(2)=-0.98415243845764617655228962221207;
xi(3)=-0.94720428399922868052421376661573;
xi(4)=-0.89006229019090447052965782577909;
xi(5)=-0.81394892761192113604544184805614;
xi(6)=-0.72048723996120215811988189639847;
xi(7)=-0.61166943828425897122621160586993;
xi(8)=-0.48981487518990234980875123568327;
xi(9)=-0.35752071013891953806095728024018;
xi(10)=-0.21760658515928504178795509346539;
xi(11)=-0.073054540010898334761088790464107;
xi(12)=0.073054540010898334761088790464107;
xi(13)=0.21760658515928504178795509346539;
xi(14)=0.35752071013891953806095728024018;
xi(15)=0.48981487518990234980875123568327;
xi(16)=0.61166943828425897122621160586993;
xi(17)=0.72048723996120215811988189639847;
xi(18)=0.81394892761192113604544184805614;
xi(19)=0.89006229019090447052965782577909;
xi(20)=0.94720428399922868052421376661573;
xi(21)=0.98415243845764617655228962221207;
xi(22)=1.0000000000000000000000000000000;
case(23)
xi(1)=-1.0000000000000000000000000000000;
xi(2)=-0.98552715587873257808146276673810;
xi(3)=-0.95175795571071020413563967985143;
xi(4)=-0.89945855804034501095016032034737;
xi(5)=-0.82965109665128588622320061929000;
xi(6)=-0.74369504117206068394516354306700;
xi(7)=-0.64326364446013620847614553360277;
xi(8)=-0.53031177113684416813011532015230;
xi(9)=-0.40703793791447482919595048821510;
xi(10)=-0.27584154894579306710687763267914;
xi(11)=-0.13927620404066839859186261298277;
xi(12)=0.0;
xi(13)=0.13927620404066839859186261298277;
xi(14)=0.27584154894579306710687763267914;
xi(15)=0.40703793791447482919595048821510;
xi(16)=0.53031177113684416813011532015230;
xi(17)=0.64326364446013620847614553360277;
xi(18)=0.74369504117206068394516354306700;
xi(19)=0.82965109665128588622320061929000;
xi(20)=0.89945855804034501095016032034737;
xi(21)=0.95175795571071020413563967985143;
xi(22)=0.98552715587873257808146276673810;
xi(23)=1.0000000000000000000000000000000;
case(24)
xi(1)=-1.0000000000000000000000000000000;
xi(2)=-0.98673055350516088355308673815447;
xi(3)=-0.95574822092988635802697713055064;
xi(4)=-0.90770567511350652199515299646621;
xi(5)=-0.84346407015487204062330503742334;
xi(6)=-0.76417048242049330778737528095229;
xi(7)=-0.67124010526412869983566485818701;
xi(8)=-0.56633135797929531218940954454228;
xi(9)=-0.45131637321432261824821849156962;
xi(10)=-0.32824761337551091203338917935961;
xi(11)=-0.19932125339083266723657253912499;
xi(12)=-0.066837993737228578113641808391677;
xi(13)=0.066837993737228578113641808391677;
xi(14)=0.19932125339083266723657253912499;
xi(15)=0.32824761337551091203338917935961;
xi(16)=0.45131637321432261824821849156962;
xi(17)=0.56633135797929531218940954454228;
xi(18)=0.67124010526412869983566485818701;
xi(19)=0.76417048242049330778737528095229;
xi(20)=0.84346407015487204062330503742334;
xi(21)=0.90770567511350652199515299646621;
xi(22)=0.95574822092988635802697713055064;
xi(23)=0.98673055350516088355308673815447;
xi(24)=1.0000000000000000000000000000000;
case(25)
xi(1)=-1.0000000000000000000000000000000;
xi(2)=-0.98778994493149370927180407087316;
xi(3)=-0.95926413825253447885984517883952;
xi(4)=-0.91498277073462257832314993373685;
xi(5)=-0.85567646583531657752381326017175;
xi(6)=-0.78231965924071678039918795222875;
xi(7)=-0.69611704881513436676036543378938;
xi(8)=-0.59848414727999326809762099718107;
xi(9)=-0.49102411481887838261895992256880;
xi(10)=-0.37550145785922723322871461226072;
xi(11)=-0.25381306416887658017988688168728;
xi(12)=-0.12795705948310697270898462509465;
xi(13)=0.0;
xi(14)=0.12795705948310697270898462509465;
xi(15)=0.25381306416887658017988688168728;
xi(16)=0.37550145785922723322871461226072;
xi(17)=0.49102411481887838261895992256880;
xi(18)=0.59848414727999326809762099718107;
xi(19)=0.69611704881513436676036543378938;
xi(20)=0.78231965924071678039918795222875;
xi(21)=0.85567646583531657752381326017175;
xi(22)=0.91498277073462257832314993373685;
xi(23)=0.95926413825253447885984517883952;
xi(24)=0.98778994493149370927180407087316;
xi(25)=1.0000000000000000000000000000000;
case(26)
xi(1)=-1.0000000000000000000000000000000;
xi(2)=-0.98872741231147565442924669037130;
xi(3)=-0.96237787476771732975066267128213;
xi(4)=-0.92143554681755734112513052174869;
xi(5)=-0.86652432395912357097445149780977;
xi(6)=-0.79847718310743743903678763857034;
xi(7)=-0.71832581636266508052540538790393;
xi(8)=-0.62728529949231688780051890716122;
xi(9)=-0.52673574202987854529824339438099;
xi(10)=-0.41820138706624678556368388233385;
xi(11)=-0.30332751285925272077404103790211;
xi(12)=-0.18385549527005490132120567188576;
xi(13)=-0.061596411781919728205487274781751;
xi(14)=0.061596411781919728205487274781751;
xi(15)=0.18385549527005490132120567188576;
xi(16)=0.30332751285925272077404103790211;
xi(17)=0.41820138706624678556368388233385;
xi(18)=0.52673574202987854529824339438099;
xi(19)=0.62728529949231688780051890716122;
xi(20)=0.71832581636266508052540538790393;
xi(21)=0.79847718310743743903678763857034;
xi(22)=0.86652432395912357097445149780977;
xi(23)=0.92143554681755734112513052174869;
xi(24)=0.96237787476771732975066267128213;
xi(25)=0.98872741231147565442924669037130;
xi(26)=1.0000000000000000000000000000000;
case(27)
xi(1)=-1.0000000000000000000000000000000;
xi(2)=-0.98956096372855062114887464005452;
xi(3)=-0.96514840245081891926426680312113;
xi(4)=-0.92718345872511577885974676407272;
xi(5)=-0.87620208621452247866654304919414;
xi(6)=-0.81292048689581228825809186475854;
xi(7)=-0.73822714984645991083942722574841;
xi(8)=-0.65317066369680951517870487518246;
xi(9)=-0.55894506094256118006266807486812;
xi(10)=-0.45687307561408240367861265083924;
xi(11)=-0.34838758198902867043708657858717;
xi(12)=-0.23501148310291813333680158051333;
xi(13)=-0.11833633389852104843870826827532;
xi(14)=0.0;
xi(15)=0.11833633389852104843870826827532;
xi(16)=0.23501148310291813333680158051333;
xi(17)=0.34838758198902867043708657858717;
xi(18)=0.45687307561408240367861265083924;
xi(19)=0.55894506094256118006266807486812;
xi(20)=0.65317066369680951517870487518246;
xi(21)=0.73822714984645991083942722574841;
xi(22)=0.81292048689581228825809186475854;
xi(23)=0.87620208621452247866654304919414;
xi(24)=0.92718345872511577885974676407272;
xi(25)=0.96514840245081891926426680312113;
xi(26)=0.98956096372855062114887464005452;
xi(27)=1.0000000000000000000000000000000;
case(28)
xi(1)=-1.0000000000000000000000000000000;
xi(2)=-0.99030540261845412551709291186460;
xi(3)=-0.96762428585713130434377087914347;
xi(4)=-0.93232516712155852492614241860899;
xi(5)=-0.88487101721130283841230633825725;
xi(6)=-0.82588097005633819604451442619660;
xi(7)=-0.75612419400556976515209572526619;
xi(8)=-0.67651012892957331771938124507339;
xi(9)=-0.58807668983717560647831207976171;
xi(10)=-0.49197675393157938042961668294448;
xi(11)=-0.38946313757636280793588415520411;
xi(12)=-0.28187226662160237195415250890933;
xi(13)=-0.17060675530800436088270267774877;
xi(14)=-0.057117121693512897626543651573747;
xi(15)=0.057117121693512897626543651573747;
xi(16)=0.17060675530800436088270267774877;
xi(17)=0.28187226662160237195415250890933;
xi(18)=0.38946313757636280793588415520411;
xi(19)=0.49197675393157938042961668294448;
xi(20)=0.58807668983717560647831207976171;
xi(21)=0.67651012892957331771938124507339;
xi(22)=0.75612419400556976515209572526619;
xi(23)=0.82588097005633819604451442619660;
xi(24)=0.88487101721130283841230633825725;
xi(25)=0.93232516712155852492614241860899;
xi(26)=0.96762428585713130434377087914347;
xi(27)=0.99030540261845412551709291186460;
xi(28)=1.0000000000000000000000000000000;
case(29)
xi(1)=-1.0000000000000000000000000000000;
xi(2)=-0.99097298826856981627874511329957;
xi(3)=-0.96984580728793629371240138788782;
xi(4)=-0.93694271852098251154203684828336;
xi(5)=-0.89266571997608827876422991405309;
xi(6)=-0.83755273628178656524207331618682;
xi(7)=-0.77227289720646874135417063413645;
xi(8)=-0.69761866135636798839397442857709;
xi(9)=-0.61449625220343294997090361742436;
xi(10)=-0.52391467437196908357881657312361;
xi(11)=-0.42697347171349438763757641096836;
xi(12)=-0.32484938284191100626001402307884;
xi(13)=-0.21878205828426101980905185232499;
xi(14)=-0.11005901339559210970778760225568;
xi(15)=0.0;
xi(16)=0.11005901339559210970778760225568;
xi(17)=0.21878205828426101980905185232499;
xi(18)=0.32484938284191100626001402307884;
xi(19)=0.42697347171349438763757641096836;
xi(20)=0.52391467437196908357881657312361;
xi(21)=0.61449625220343294997090361742436;
xi(22)=0.69761866135636798839397442857709;
xi(23)=0.77227289720646874135417063413645;
xi(24)=0.83755273628178656524207331618682;
xi(25)=0.89266571997608827876422991405309;
xi(26)=0.93694271852098251154203684828336;
xi(27)=0.96984580728793629371240138788782;
xi(28)=0.99097298826856981627874511329957;
xi(29)=1.0000000000000000000000000000000;
case(30)
xi(1)=-1.0000000000000000000000000000000;
xi(2)=-0.99157394284050029333881822590865;
xi(3)=-0.97184660316626924167659295220919;
xi(4)=-0.94110478095105708230718954941552;
xi(5)=-0.89969921819927685955334282314945;
xi(6)=-0.84809948718019810955142383218307;
xi(7)=-0.78689035723754708044995952059240;
xi(8)=-0.71676539863708513163379859612990;
xi(9)=-0.63851917580755840737109348270974;
xi(10)=-0.55303826009505285238455023004066;
xi(11)=-0.46129119016824068522655722110715;
xi(12)=-0.36431750042244899775598510555179;
xi(13)=-0.26321594371957379126709946900007;
xi(14)=-0.15913204262585046782494857709885;
xi(15)=-0.053245110485486669363007834451354;
xi(16)=0.053245110485486669363007834451354;
xi(17)=0.15913204262585046782494857709885;
xi(18)=0.26321594371957379126709946900007;
xi(19)=0.36431750042244899775598510555179;
xi(20)=0.46129119016824068522655722110715;
xi(21)=0.55303826009505285238455023004066;
xi(22)=0.63851917580755840737109348270974;
xi(23)=0.71676539863708513163379859612990;
xi(24)=0.78689035723754708044995952059240;
xi(25)=0.84809948718019810955142383218307;
xi(26)=0.89969921819927685955334282314945;
xi(27)=0.94110478095105708230718954941552;
xi(28)=0.97184660316626924167659295220919;
xi(29)=0.99157394284050029333881822590865;
xi(30)=1.0000000000000000000000000000000;
case(31)
xi(1)=-1.0000000000000000000000000000000;
xi(2)=-0.99211684434648108858932461564871;
xi(3)=-0.97365493581573647482606954847523;
xi(4)=-0.94486917020803922980640553764986;
xi(5)=-0.90606695144126981314451152324746;
xi(6)=-0.85765999529745561036700538239647;
xi(7)=-0.80016154319246298927233133398318;
xi(8)=-0.73418113630907524213414083483764;
xi(9)=-0.66041820261152552333704127459677;
xi(10)=-0.57965465720800191300329686219773;
xi(11)=-0.49274661909883240570380023307021;
xi(12)=-0.40061533828056180630266441939609;
xi(13)=-0.30423743127287262095945514762626;
xi(14)=-0.20463452924752450114109358572771;
xi(15)=-0.10286244876068227708911472330805;
xi(16)=0.0;
xi(17)=0.10286244876068227708911472330805;
xi(18)=0.20463452924752450114109358572771;
xi(19)=0.30423743127287262095945514762626;
xi(20)=0.40061533828056180630266441939609;
xi(21)=0.49274661909883240570380023307021;
xi(22)=0.57965465720800191300329686219773;
xi(23)=0.66041820261152552333704127459677;
xi(24)=0.73418113630907524213414083483764;
xi(25)=0.80016154319246298927233133398318;
xi(26)=0.85765999529745561036700538239647;
xi(27)=0.90606695144126981314451152324746;
xi(28)=0.94486917020803922980640553764986;
xi(29)=0.97365493581573647482606954847523;
xi(30)=0.99211684434648108858932461564871;
xi(31)=1.0000000000000000000000000000000;
case(32)
xi(1)=-1.0000000000000000000000000000000;
xi(2)=-0.99260893397276135937154131527745;
xi(3)=-0.97529469048270922806245595038815;
xi(4)=-0.94828483841723237808332824158843;
xi(5)=-0.91184993906373190407453145237407;
xi(6)=-0.86635247601267551983081476942454;
xi(7)=-0.81224473177744234454682178807258;
xi(8)=-0.75006449393667479771727596467830;
xi(9)=-0.68042975561555081594239575769885;
xi(10)=-0.60403258714842112613717291261189;
xi(11)=-0.52163226288156529060659257339202;
xi(12)=-0.43404771720184693960331935700639;
xi(13)=-0.34214940653888148625381222069872;
xi(14)=-0.24685065885020530441624143115887;
xi(15)=-0.14909859681364749491438187617004;
xi(16)=-0.049864725046593252306297744728426;
xi(17)=0.049864725046593252306297744728426;
xi(18)=0.14909859681364749491438187617004;
xi(19)=0.24685065885020530441624143115887;
xi(20)=0.34214940653888148625381222069872;
xi(21)=0.43404771720184693960331935700639;
xi(22)=0.52163226288156529060659257339202;
xi(23)=0.60403258714842112613717291261189;
xi(24)=0.68042975561555081594239575769885;
xi(25)=0.75006449393667479771727596467830;
xi(26)=0.81224473177744234454682178807258;
xi(27)=0.86635247601267551983081476942454;
xi(28)=0.91184993906373190407453145237407;
xi(29)=0.94828483841723237808332824158843;
xi(30)=0.97529469048270922806245595038815;
xi(31)=0.99260893397276135937154131527745;
xi(32)=1.0000000000000000000000000000000; case(33)
xi(1)=-1.0000000000000000000000000000000;
xi(2)=-0.99305635843365834366733608675624;
xi(3)=-0.97678616331690630148562760838132;
xi(4)=-0.95139345139699574337552441691356;
xi(5)=-0.91711730345094124082567906975682;
xi(6)=-0.87427810075056222064656062933034;
xi(7)=-0.82327592300406746955907131326671;
xi(8)=-0.76458700179352862800665848042309;
xi(9)=-0.69875931661816259569967724754442;
xi(10)=-0.62640749128126825726383134241071;
xi(11)=-0.54820705991911162311040471176063;
xi(12)=-0.46488816163210675597300372098526;
xi(13)=-0.37722872425339363504193466071800;
xi(14)=-0.2860472014876740416053445029630;
xi(15)=-0.19219493146747722574130699544116;
xi(16)=-0.096548188176107006316957788091625;
xi(17)=0.0;
xi(18)=0.096548188176107006316957788091625;
xi(19)=0.19219493146747722574130699544116;
xi(20)=0.28604720148767404160534450296304;
xi(21)=0.37722872425339363504193466071800;
xi(22)=0.46488816163210675597300372098526;
xi(23)=0.54820705991911162311040471176063;
xi(24)=0.62640749128126825726383134241071;
xi(25)=0.69875931661816259569967724754442;
xi(26)=0.76458700179352862800665848042309;
xi(27)=0.82327592300406746955907131326671;
xi(28)=0.87427810075056222064656062933034;
xi(29)=0.91711730345094124082567906975682;
xi(30)=0.95139345139699574337552441691356;
xi(31)=0.97678616331690630148562760838132;
xi(32)=0.99305635843365834366733608675624;
xi(33)=1.0000000000000000000000000000000;
otherwise
error('LOBATTO_POINTS - Fatal error! Illegal value of n.');
end
return
end
function w = lobatto_weights(n)
w = ones(1,n);
switch n
case(2)
w(1)=1.0000000000000000000000000000000;
w(2)=1.0000000000000000000000000000000;
case(3)
w(1)=0.33333333333333333333333333333333;
w(2)=1.3333333333333333333333333333333;
w(3)=0.33333333333333333333333333333333;
case(4)
w(1)=0.16666666666666666666666666666667;
w(2)=0.83333333333333333333333333333333;
w(3)=0.83333333333333333333333333333333;
w(4)=0.16666666666666666666666666666667;
case(5)
w(1)=0.10000000000000000000000000000000;
w(2)=0.54444444444444444444444444444444;
w(3)=0.71111111111111111111111111111111;
w(4)=0.54444444444444444444444444444444;
w(5)=0.10000000000000000000000000000000;
case(6)
w(1)=0.066666666666666666666666666666667;
w(2)=0.37847495629784698031661280821202;
w(3)=0.55485837703548635301672052512131;
w(4)=0.55485837703548635301672052512131;
w(5)=0.37847495629784698031661280821202;
w(6)=0.066666666666666666666666666666667;
case(7)
w(1)=0.047619047619047619047619047619048;
w(2)=0.27682604736156594801070040629007;
w(3)=0.43174538120986262341787102228136;
w(4)=0.48761904761904761904761904761905;
w(5)=0.43174538120986262341787102228136;
w(6)=0.27682604736156594801070040629007;
w(7)=0.047619047619047619047619047619048;
case(8)
w(1)=0.035714285714285714285714285714286;
w(2)=0.21070422714350603938299206577576;
w(3)=0.34112269248350436476424067710775;
w(4)=0.41245879465870388156705297140221;
w(5)=0.41245879465870388156705297140221;
w(6)=0.34112269248350436476424067710775;
w(7)=0.21070422714350603938299206577576;
w(8)=0.035714285714285714285714285714286;
case(9)
w(1)=0.027777777777777777777777777777778;
w(2)=0.16549536156080552504633972002921;
w(3)=0.27453871250016173528070561857937;
w(4)=0.34642851097304634511513153213972;
w(5)=0.37151927437641723356009070294785;
w(6)=0.34642851097304634511513153213972;
w(7)=0.27453871250016173528070561857937;
w(8)=0.16549536156080552504633972002921;
w(9)=0.027777777777777777777777777777778;
case(10)
w(1)=0.022222222222222222222222222222222;
w(2)=0.13330599085107011112622717075539;
w(3)=0.22488934206312645211945782173105;
w(4)=0.29204268367968375787558225737444;
w(5)=0.32753976118389745665651052791689;
w(6)=0.32753976118389745665651052791689;
w(7)=0.29204268367968375787558225737444;
w(8)=0.22488934206312645211945782173105;
w(9)=0.13330599085107011112622717075539
w(10)=0.022222222222222222222222222222222;
case(11)
w(1)=0.018181818181818181818181818181818;
w(2)=0.10961227326699486446140344958035;
w(3)=0.18716988178030520410814152189943;
w(4)=0.24804810426402831404008486642187;
w(5)=0.28687912477900808867922240333154;
w(6)=0.30021759545569069378593188116998;
w(7)=0.28687912477900808867922240333154;
w(8)=0.24804810426402831404008486642187;
w(9)=0.18716988178030520410814152189943;
w(10)=0.10961227326699486446140344958035;
w(11)=0.018181818181818181818181818181818;
case(12)
w(1)=0.015151515151515151515151515151515;
w(2)=0.091684517413196130668342594134079;
w(3)=0.15797470556437011516467106270034;
w(4)=0.21250841776102114535830207736687;
w(5)=0.25127560319920128029324441214760;
w(6)=0.27140524091069617700028833849960;
w(7)=0.27140524091069617700028833849960;
w(8)=0.25127560319920128029324441214760;
w(9)=0.21250841776102114535830207736687;
w(10)=0.15797470556437011516467106270034;
w(11)=0.091684517413196130668342594134079;
w(12)=0.015151515151515151515151515151515;
case(13)
w(1)=0.012820512820512820512820512820513;
w(2)=0.077801686746818927793588988333134;
w(3)=0.13498192668960834911991476258937;
w(4)=0.18364686520355009200749425874681;
w(5)=0.22076779356611008608553400837940;
w(6)=0.24401579030667635645857814836016;
w(7)=0.25193084933344673604413864154124;
w(8)=0.24401579030667635645857814836016;
w(9)=0.22076779356611008608553400837940;
w(10)=0.18364686520355009200749425874681;
w(11)=0.13498192668960834911991476258937;
w(12)=0.077801686746818927793588988333134;
w(13)=0.012820512820512820512820512820513;
case(14)
w(1)=0.010989010989010989010989010989011;
w(2)=0.066837284497681284634070660746053;
w(3)=0.11658665589871165154099667065465;
w(4)=0.16002185176295214241282099798759;
w(5)=0.19482614937341611864033177837588;
w(6)=0.21912625300977075487116252395417;
w(7)=0.23161279446845705888962835729264;
w(8)=0.23161279446845705888962835729264;
w(9)=0.21912625300977075487116252395417;
w(10)=0.19482614937341611864033177837588;
w(11)=0.16002185176295214241282099798759;
w(12)=0.11658665589871165154099667065465;
w(13)=0.066837284497681284634070660746053;
w(14)=0.010989010989010989010989010989011;
case(15)
w(1)=0.0095238095238095238095238095238095;
w(2)=0.058029893028601249096880584025282;
w(3)=0.10166007032571806760366617078880;
w(4)=0.14051169980242810946044680564367;
w(5)=0.17278964725360094905207709940835;
w(6)=0.19698723596461335609250034650741;
w(7)=0.21197358592682092012743007697722;
w(8)=0.21704811634881564951495021425091;
w(9)=0.21197358592682092012743007697722;
w(10)=0.19698723596461335609250034650741;
w(11)=0.17278964725360094905207709940835;
w(12)=0.14051169980242810946044680564367;
w(13)=0.10166007032571806760366617078880;
w(14)=0.058029893028601249096880584025282;
w(15)=0.0095238095238095238095238095238095;
case(16)
w(1)=0.0083333333333333333333333333333333;
w(2)=0.050850361005919905403244919565455;
w(3)=0.089393697325930800991052080166084;
w(4)=0.12425538213251409834953633265731;
w(5)=0.15402698080716428081564494048499;
w(6)=0.17749191339170412530107566952836;
w(7)=0.19369002382520358431691359885352;
w(8)=0.20195830817822987148919912541094;
w(9)=0.20195830817822987148919912541094;
w(10)=0.19369002382520358431691359885352;
w(11)=0.17749191339170412530107566952836;
w(12)=0.15402698080716428081564494048499;
w(13)=0.12425538213251409834953633265731;
w(14)=0.089393697325930800991052080166084;
w(15)=0.050850361005919905403244919565455;
w(16)=0.0083333333333333333333333333333333;
case(17)
w(1)=0.0073529411764705882352941176470588;
w(2)=0.044921940543254209647400954623212;
w(3)=0.079198270503687119190264429952835;
w(4)=0.11059290900702816137577270522008;
w(5)=0.13798774620192655905620157495403;
w(6)=0.16039466199762153951632836586475;
w(7)=0.17700425351565787043694574536329;
w(8)=0.18721633967761923589208848286062;
w(9)=0.19066187475346943329940724702825;
w(10)=0.18721633967761923589208848286062;
w(11)=0.17700425351565787043694574536329;
w(12)=0.16039466199762153951632836586475;
w(13)=0.13798774620192655905620157495403;
w(14)=0.11059290900702816137577270522008;
w(15)=0.079198270503687119190264429952835;
w(16)=0.044921940543254209647400954623212;
w(17)=0.0073529411764705882352941176470588;
case(18)
w(1)=0.0065359477124183006535947712418301;
w(2)=0.039970628810914066137599176410101;
w(3)=0.070637166885633664999222960167786;
w(4)=0.099016271717502802394423605318672;
w(5)=0.12421053313296710026339635889675;
w(6)=0.14541196157380226798300321049443;
w(7)=0.16193951723760248926432670670023;
w(8)=0.17326210948945622601061440382668;
w(9)=0.17901586343970308229381880694353;
w(10)=0.17901586343970308229381880694353;
w(11)=0.17326210948945622601061440382668;
w(12)=0.16193951723760248926432670670023;
w(13)=0.14541196157380226798300321049443;
w(14)=0.12421053313296710026339635889675;
w(15)=0.099016271717502802394423605318672;
w(16)=0.070637166885633664999222960167786;
w(17)=0.039970628810914066137599176410101;
w(18)=0.0065359477124183006535947712418301;
case(19)
w(1)=0.0058479532163742690058479532163743
w(2)=0.035793365186176477115425569035122
w(3)=0.063381891762629736851695690418317
w(4)=0.089131757099207084448008790556153
w(5)=0.11231534147730504407091001546378
w(6)=0.13226728044875077692604673390973
w(7)=0.14841394259593888500968064366841
w(8)=0.16029092404406124197991096818359
w(9)=0.16755658452714286727013727774026
w(10)=0.17000191928482723464467271561652
w(11)=0.16755658452714286727013727774026
w(12)=0.16029092404406124197991096818359
w(13)=0.14841394259593888500968064366841
w(14)=0.13226728044875077692604673390973
w(15)=0.11231534147730504407091001546378
w(16)=0.089131757099207084448008790556153
w(17)=0.063381891762629736851695690418317
w(18)=0.035793365186176477115425569035122
w(19)=0.0058479532163742690058479532163743
case(20)
w(1)=0.0052631578947368421052631578947368
w(2)=0.032237123188488941491605028117294
w(3)=0.057181802127566826004753627173243
w(4)=0.080631763996119603144776846113721
w(5)=0.10199149969945081568378120573289
w(6)=0.12070922762867472509942970500239
w(7)=0.13630048235872418448978079298903
w(8)=0.14836155407091682581471301373397
w(9)=0.15658010264747548715816989679364
w(10)=0.16074328638784574900772672644908
w(11)=0.16074328638784574900772672644908
w(12)=0.15658010264747548715816989679364
w(13)=0.14836155407091682581471301373397
w(14)=0.13630048235872418448978079298903
w(15)=0.12070922762867472509942970500239
w(16)=0.10199149969945081568378120573289
w(17)=0.080631763996119603144776846113721
w(18)=0.057181802127566826004753627173243
w(19)=0.032237123188488941491605028117294
w(20)=0.0052631578947368421052631578947368
case(21)
w(1)=0.0047619047619047619047619047619048
w(2)=0.029184840098505458609458543613171
w(3)=0.051843169000849625072722971852830
w(4)=0.073273918185074144252547861041894
w(5)=0.092985467957886065301137664149214
w(6)=0.11051708321912333526700048678439
w(7)=0.12545812119086894801515753570800
w(8)=0.13745846286004134358089961741515
w(9)=0.14623686244797745926727053063439
w(10)=0.15158757511168138445325068150529
w(11)=0.15338519033217494855158440506754
w(12)=0.15158757511168138445325068150529
w(13)=0.14623686244797745926727053063439
w(14)=0.13745846286004134358089961741515
w(15)=0.12545812119086894801515753570800
w(16)=0.11051708321912333526700048678439
w(17)=0.092985467957886065301137664149214
w(18)=0.073273918185074144252547861041894
w(19)=0.051843169000849625072722971852830
w(20)=0.029184840098505458609458543613171
w(21)=0.0047619047619047619047619047619048
case(22)
w(1)=0.0043290043290043290043290043290043
w(2)=0.026545747682501757911627904520543
w(3)=0.047214465293740752123775734864792
w(4)=0.066865605864553076012404194157097
w(5)=0.085090060391838447815711236095748
w(6)=0.10150057480164767437243730374960
w(7)=0.11574764465393906659003636772146
w(8)=0.12752769665343027553084445930883
w(9)=0.13658968861374142668617736220617
w(10)=0.14274049227136140033623599356679
w(11)=0.14584901944424179361642043947997
w(12)=0.14584901944424179361642043947997
w(13)=0.14274049227136140033623599356679
w(14)=0.13658968861374142668617736220617
w(15)=0.12752769665343027553084445930883
w(16)=0.11574764465393906659003636772146
w(17)=0.10150057480164767437243730374960
w(18)=0.085090060391838447815711236095748
w(19)=0.066865605864553076012404194157097
w(20)=0.047214465293740752123775734864792
w(21)=0.026545747682501757911627904520543
w(22)=0.0043290043290043290043290043290043
case(23)
w(1)=0.0039525691699604743083003952569170
w(2)=0.024248600771531736517399658937097
w(3)=0.043175871170241834748876465612042
w(4)=0.061252477129554206381382847440355
w(5)=0.078135449475569989741934255347965
w(6)=0.093497246163512341833500706906697
w(7)=0.10703910172433651153518362791547
w(8)=0.11849751066274913130212600472426
w(9)=0.12764947470175887663614855305567
w(10)=0.13431687263860381990156489770071
w(11)=0.13836993638580739452350273386294
w(12)=0.13972978001274736514015970647975
w(13)=0.13836993638580739452350273386294
w(14)=0.13431687263860381990156489770071
w(15)=0.12764947470175887663614855305567
w(16)=0.11849751066274913130212600472426
w(17)=0.10703910172433651153518362791547
w(18)=0.093497246163512341833500706906697
w(19)=0.078135449475569989741934255347965
w(20)=0.061252477129554206381382847440355
w(21)=0.043175871170241834748876465612042
w(22)=0.024248600771531736517399658937097
w(23)=0.0039525691699604743083003952569170
case(24)
w(1)=0.0036231884057971014492753623188406
w(2)=0.022236853464711208992960434817299
w(3)=0.039631681333467809469262163363472
w(4)=0.056309848724646199020948339678118
w(5)=0.071981862055293982215639060668680
w(6)=0.086369029967929068216562819306583
w(7)=0.099214827684083587414149770995712
w(8)=0.11029008689296860411042412302760
w(9)=0.11939719370249131903229388835072
w(10)=0.12637364202802080012724809069137
w(11)=0.13109494187360394235445022581653
w(12)=0.13347684386698637759678572096508
w(13)=0.13347684386698637759678572096508
w(14)=0.13109494187360394235445022581653
w(15)=0.12637364202802080012724809069137
w(16)=0.11939719370249131903229388835072
w(17)=0.11029008689296860411042412302760
w(18)=0.099214827684083587414149770995712
w(19)=0.086369029967929068216562819306583
w(20)=0.071981862055293982215639060668680
w(21)=0.056309848724646199020948339678118
w(22)=0.039631681333467809469262163363472
w(23)=0.022236853464711208992960434817299
w(24)=0.0036231884057971014492753623188406
case(25)
w(1)=0.0033333333333333333333333333333333
w(2)=0.020465168932974385308542471819163
w(3)=0.036504738794271372032382988755110
w(4)=0.051936228368491474643333889713489
w(5)=0.066513728675312784693869993313599
w(6)=0.079998774836292981801626436138836
w(7)=0.092170139910620421912689622714049
w(8)=0.10282803034795783082750364490713
w(9)=0.11179746626832088815624423243104
w(10)=0.11893117940681182540944424463419
w(11)=0.12411203893795029069521531243930
w(12)=0.12725497753833144701711441673876
w(13)=0.12830838929866192833739882612401
w(14)=0.12725497753833144701711441673876
w(15)=0.12411203893795029069521531243930
w(16)=0.11893117940681182540944424463419
w(17)=0.11179746626832088815624423243104
w(18)=0.10282803034795783082750364490713
w(19)=0.092170139910620421912689622714049
w(20)=0.079998774836292981801626436138836
w(21)=0.066513728675312784693869993313599
w(22)=0.051936228368491474643333889713489
w(23)=0.036504738794271372032382988755110
w(24)=0.020465168932974385308542471819163
w(25)=0.0033333333333333333333333333333333
case(26)
w(1)=0.0030769230769230769230769230769231
w(2)=0.018896858024263465581352047676267
w(3)=0.033732303685955999377522748582756
w(4)=0.048048399081180627315979312780662
w(5)=0.061635025142547402782180054104279
w(6)=0.074287050122291137316084082986249
w(7)=0.085812863980004362187481252359068
w(8)=0.096037802353901310803697396370862
w(9)=0.10480688623073705298990583689136
w(10)=0.11198719411986033530021087987989
w(11)=0.11746988409380900707669562029648
w(12)=0.12117184628844334604428669570423
w(13)=0.12303696380008287630152714929097
w(14)=0.12303696380008287630152714929097
w(15)=0.12117184628844334604428669570423
w(16)=0.11746988409380900707669562029648
w(17)=0.11198719411986033530021087987989
w(18)=0.10480688623073705298990583689136
w(19)=0.096037802353901310803697396370862
w(20)=0.085812863980004362187481252359068
w(21)=0.074287050122291137316084082986249
w(22)=0.061635025142547402782180054104279
w(23)=0.048048399081180627315979312780662
w(24)=0.033732303685955999377522748582756
w(25)=0.018896858024263465581352047676267
w(26)=0.0030769230769230769230769230769231
case(27)
w(1)=0.0028490028490028490028490028490028
w(2)=0.017501974876065579019369734362719
w(3)=0.031262951735202384324784990064795
w(4)=0.044577657933061698744074742241386
w(5)=0.057265569680162731739087228253941
w(6)=0.069149342360043276280634735030127
w(7)=0.080062321970538458168238070483950
w(8)=0.089851365259290559972136301946063
w(9)=0.098379074585952763179308931269989
w(10)=0.10552574782125301121838198947990
w(11)=0.11119106525743703292826195934215
w(12)=0.11529550025465198281375333337473
w(13)=0.11778143658595615906649907944836
w(14)=0.11861397766276302708523980370575
w(15)=0.11778143658595615906649907944836
w(16)=0.11529550025465198281375333337473
w(17)=0.11119106525743703292826195934215
w(18)=0.10552574782125301121838198947990
w(19)=0.098379074585952763179308931269989
w(20)=0.089851365259290559972136301946063
w(21)=0.080062321970538458168238070483950
w(22)=0.069149342360043276280634735030127
w(23)=0.057265569680162731739087228253941
w(24)=0.044577657933061698744074742241386
w(25)=0.031262951735202384324784990064795
w(26)=0.017501974876065579019369734362719
w(27)=0.0028490028490028490028490028490028
case(28)
w(1)=0.0026455026455026455026455026455026
w(2)=0.016255883957504218198990152897110
w(3)=0.029054220677979144706852779823582
w(4)=0.041466915243006721037985848778627
w(5)=0.053338077047327488400325513471500
w(6)=0.064513658080354538390245693462478
w(7)=0.074848123509707702780567939143162
w(8)=0.084206795121510325739491187797920
w(9)=0.092467685997712175722035394699386
w(10)=0.099523110412495672875945354750915
w(11)=0.10528109376105589448237755179491
w(12)=0.10966657379597662508858586379121
w(13)=0.11262238007723901254063642411366
w(14)=0.11410997967262783453331479283004
w(15)=0.11410997967262783453331479283004
w(16)=0.11262238007723901254063642411366
w(17)=0.10966657379597662508858586379121
w(18)=0.10528109376105589448237755179491
w(19)=0.099523110412495672875945354750915
w(20)=0.092467685997712175722035394699386
w(21)=0.084206795121510325739491187797920
w(22)=0.074848123509707702780567939143162
w(23)=0.064513658080354538390245693462478
w(24)=0.053338077047327488400325513471500
w(25)=0.041466915243006721037985848778627
w(26)=0.029054220677979144706852779823582
w(27)=0.016255883957504218198990152897110
w(28)=0.0026455026455026455026455026455026
case(29)
w(1)=0.0024630541871921182266009852216749
w(2)=0.015138169859967596790857926606856
w(3)=0.027070806296824827540851425445026
w(4)=0.038668439979712979747832574267124
w(5)=0.049795809093237560922668341175478
w(6)=0.060318503828522735028064463517634
w(7)=0.070108938000597999185266858424014
w(8)=0.079048313027885601784170232175109
w(9)=0.087028133441135605637626155215633
w(10)=0.093951542114796312497836506987570
w(11)=0.099734501637149860653569177534068
w(12)=0.10430681646367653742270711219124
w(13)=0.10761298583356853937358685404661
w(14)=0.10961287784589042013213241075660
w(15)=0.11028221677968261011245795287074
w(16)=0.10961287784589042013213241075660
w(17)=0.10761298583356853937358685404661
w(18)=0.10430681646367653742270711219124
w(19)=0.099734501637149860653569177534068
w(20)=0.093951542114796312497836506987570
w(21)=0.087028133441135605637626155215633
w(22)=0.079048313027885601784170232175109
w(23)=0.070108938000597999185266858424014
w(24)=0.060318503828522735028064463517634
w(25)=0.049795809093237560922668341175478
w(26)=0.038668439979712979747832574267124
w(27)=0.027070806296824827540851425445026
w(28)=0.015138169859967596790857926606856
w(29)=0.0024630541871921182266009852216749
case(30)
w(1)=0.0022988505747126436781609195402299
w(2)=0.014131799327905387640732168668208
w(3)=0.025283166740551402204268253850042
w(4)=0.036142094199408535314732683123360
w(5)=0.046590694533142927401880491446022
w(6)=0.056511197923080383302193710472303
w(7)=0.065791336397790054944101374593132
w(8)=0.074326003324718253834067642000831
w(9)=0.082018512833406914799721545527160
w(10)=0.088781712319765210167256434455921
w(11)=0.094538975193860891781132715364036
w(12)=0.099225071004299830657675964954309
w(13)=0.10278690530723498947028080833538
w(14)=0.10518412159645464985621298085914
w(15)=0.10638955872366792494758230680993
w(16)=0.10638955872366792494758230680993
w(17)=0.10518412159645464985621298085914
w(18)=0.10278690530723498947028080833538
w(19)=0.099225071004299830657675964954309
w(20)=0.094538975193860891781132715364036
w(21)=0.088781712319765210167256434455921
w(22)=0.082018512833406914799721545527160
w(23)=0.074326003324718253834067642000831
w(24)=0.065791336397790054944101374593132
w(25)=0.056511197923080383302193710472303
w(26)=0.046590694533142927401880491446022
w(27)=0.036142094199408535314732683123360
w(28)=0.025283166740551402204268253850042
w(29)=0.014131799327905387640732168668208
w(30)=0.0022988505747126436781609195402299
case(31)
w(1)=0.0021505376344086021505376344086022
w(2)=0.013222471025464670302635629869835
w(3)=0.023666433230270316730045782709426
w(4)=0.033853940405224057629810846931873
w(5)=0.043681818160066912829747288083268
w(6)=0.053046465493448782774280860067640
w(7)=0.061848741290454623390748584002232
w(8)=0.069995377594100570450215353373140
w(9)=0.077400032341475618634283878417423
w(10)=0.083984220517529730665258108226376
w(11)=0.089678151045260822248142932836209
w(12)=0.094421468377857957432190489816534
w(13)=0.098163893013712757119741879401421
w(14)=0.10086575479865051339556259610285
w(15)=0.10249841359547039767590565499043
w(16)=0.10304456295320733314178496152549
w(17)=0.10249841359547039767590565499043
w(18)=0.10086575479865051339556259610285
w(19)=0.098163893013712757119741879401421
w(20)=0.094421468377857957432190489816534
w(21)=0.089678151045260822248142932836209
w(22)=0.083984220517529730665258108226376
w(23)=0.077400032341475618634283878417423
w(24)=0.069995377594100570450215353373140
w(25)=0.061848741290454623390748584002232
w(26)=0.053046465493448782774280860067640
w(27)=0.043681818160066912829747288083268
w(28)=0.033853940405224057629810846931873
w(29)=0.023666433230270316730045782709426
w(30)=0.013222471025464670302635629869835
w(31)=0.0021505376344086021505376344086022
case(32)
w(1)=0.0020161290322580645161290322580645
w(2)=0.012398106501373843788620349229117
w(3)=0.022199552889291964623832171161080
w(4)=0.031775135410915465781562278917906
w(5)=0.041034201586062723330403841719016
w(6)=0.049885271336221207011960153724443
w(7)=0.058240497248055869550798929919204
w(8)=0.066016877257154543932436331683494
w(9)=0.073137139602679032640370983569209
w(10)=0.079530525692106252292356728883049
w(11)=0.085133497949668230527527658506617
w(12)=0.089890372957357833072124789584664
w(13)=0.093753875546813813565908354145794
w(14)=0.096685608948002600560378147706395
w(15)=0.098656436540761777170651164243408
w(16)=0.099646771501276777634939084748541
w(17)=0.099646771501276777634939084748541
w(18)=0.098656436540761777170651164243408
w(19)=0.096685608948002600560378147706395
w(20)=0.093753875546813813565908354145794
w(21)=0.089890372957357833072124789584664
w(22)=0.085133497949668230527527658506617
w(23)=0.079530525692106252292356728883049
w(24)=0.073137139602679032640370983569209
w(25)=0.066016877257154543932436331683494
w(26)=0.058240497248055869550798929919204
w(27)=0.049885271336221207011960153724443
w(28)=0.041034201586062723330403841719016
w(29)=0.031775135410915465781562278917906
w(30)=0.022199552889291964623832171161080
w(31)=0.012398106501373843788620349229117
w(32)=0.0020161290322580645161290322580645
case(33)
w(1)=0.0018939393939393939393939393939394
w(2)=0.011648448392267734651222179870285
w(3)=0.020864609017603360095811664182613
w(4)=0.029881045916746477519971156454995
w(5)=0.038617814771813967563858875675217
w(6)=0.046993850461024170547973112041225
w(7)=0.054931059442626967951698841306657
w(8)=0.062355367852465305441060959167732
w(9)=0.069197469494016147559969796800326
w(10)=0.075393486923973828507129918402909
w(11)=0.080885572193455092179707696772196
w(12)=0.085622448531813132550191963745865
w(13)=0.089559889747077400661148462725561
w(14)=0.092661133442241463529935120285105
w(15)=0.094897224394591815824355414662604
w(16)=0.096247284972985461996554494277935
w(17)=0.096698710102716558960032808469672
w(18)=0.096247284972985461996554494277935
w(19)=0.094897224394591815824355414662604
w(20)=0.092661133442241463529935120285105
w(21)=0.089559889747077400661148462725561
w(22)=0.085622448531813132550191963745865
w(23)=0.080885572193455092179707696772196
w(24)=0.075393486923973828507129918402909
w(25)=0.069197469494016147559969796800326
w(26)=0.062355367852465305441060959167732
w(27)=0.054931059442626967951698841306657
w(28)=0.046993850461024170547973112041225
w(29)=0.038617814771813967563858875675217
w(30)=0.029881045916746477519971156454995
w(31)=0.020864609017603360095811664182613
w(32)=0.011648448392267734651222179870285
w(33)=0.0018939393939393939393939393939394
otherwise
error('LOBATTO_WEIGHTS - Fatal error! Illegal value of n.');
end
return
end
- lagrange.m给出各个拉格朗日节点基函数在所给点处的函数值及导数值：
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% function [phi dphi] = lagrange(xi,n)
%
% Purpose
% =======
% Compute the Lagrange shape functions and its derivatives for xi \in [-1,1]
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%计算拉格朗日节点基函数在xi处的值及其导数值
function [phi, dphi ] = lagrange(xi,n)
% n = number of shape functions (n = order+1)
xip = lobatto_points(n);%原来的坐标已经校正了，还需要再产生拉巴托点么？
xivec = xi*ones(1,n);%xi复制成一个vector
vec = xivec - xip;%计算拉巴托点到给定点的偏移
phi = zeros(1,n);
dphi = zeros(1,n);
denom = zeros(1,n);%初始化
for i = 1:n
num = 1;
den = 1.;
for j = 1:n
if (j ~= i)
num = num*vec(j);%其他点的偏移的乘积
den = den*(xip(i) - xip(j));
end
end
denom(i) = den;
phi(i) = num/den;
end
for i = 1:n
avec = vec; avec(i) = 1.;%第i个偏移置位为1
num = 0.;
for j = 1:n
if (j ~= i)
bvec = avec; bvec(j) = 1.;
num = num + prod(bvec);%i不等于j，i对j时，将偏移向量的除i和j处外的相乘
end
end
dphi(i) = num/denom(i);
end
return
end
-setupdofs.m设置自由度，返回自由点的内部编号集以及自由点对应整体的一个编号。
function [nodaldofs,dofinfo,numdofs] = setupdofs(options,coord,connect)
numel = options.numel;
numnod = options.numnod;
isboundarynode = zeros(1,numnod);
isboundarynode(1) = 1;
isboundarynode(numnod) = 1; % Dirichlet BCs on boundary nodes
nodaldofs(numnod).list = [];%初始化一个13*1的结构体，字段名称为list
numdofs = 0;
for n = 1:numnod
if (~isboundarynode(n))
numdofs = numdofs + 1; % classical DOF
nodaldofs(n).list = [numdofs];%自由点编号
dofinfo(numdofs) = n; % node number for dof equal to numdofs
%自由点对应的整体节点编号
end
end
return
end
- printelementnumber.m根据需求打印进度。
function printelementnumber(e,numel)
if (numel < 20) nbe = 1;
elseif (numel < 100) nbe = 10;
elseif (numel < 200) nbe = 20;
elseif (numel < 300) nbe = 30;
elseif (numel < 600) nbe = 50;
elseif (numel < 1500) nbe = 100;
elseif (numel < 4000) nbe = 400;
else nbe = 1000; end
if (mod(e,nbe) == 0)
fprintf('e = %g\n',e);
end
编写该程序需要注意的是，要给单元，节点，自由节点等编号号码，特别是自由节点的存放，可以用结构体的方式，使用空集占位。这样由连接矩阵和自由节点结构体，很快就能提取出每个单元对应的自由节点。
对于强制（本质）边界条件，可以在最后刚度矩阵上处理，也可以在前面不将边界点视为自由点而产生自由度。殊途同归。
扫一扫
1159
被折叠的 条评论
为什么被折叠?
到【灌水乐园】发言
