Code covered by the BSD License  

Highlights from
FELICITY

image thumbnail

FELICITY

by

 

21 Apr 2011 (Updated )

Finite ELement Implementation and Computational Interface Tool for You

Quad_On_Triangle(Order)
function [Points, Weights] = Quad_On_Triangle(Order)
%Quad_On_Triangle
%
%   This routine gives the quadrature rule for integrating on the unit `logical' triangle,
%   i.e. the triangle whose vertices are (0,0), (1,0), (0,1) in that order.  Be careful
%   with the input argument (see below).
%
%   [Points, Weights] = Quad_On_Triangle(Order);
%
%   OUTPUTS
%   -------
%   Points:
%       An N x 2 array, where N is the number of quadrature points.  Each row gives the
%       (x,y) coordinates of the point.  Note: N := Order.
%
%   Weights:
%       An N x 1 array, where N is the number of quadrature points.  Each row gives the
%       weight assigned to the point.  Note: N := Order.  WARNING: the quad weights are
%       scaled such that they sum to (1/2); this is the area of the reference triangle.
%
%   INPUTS
%   ------
%   Order:
%       An positive integer indicating the desired order of accuracy for the integration
%       method, i.e. (h^Order) accuracy.  Note: not all orders are available;  only these:
%       Order = [1, 3, 4, 6, 7, 9, 12, 13, 19, 28, 37].

% Copyright (c) 02-08-2008,  Shawn W. Walker

% make sure it is an integer
Order = round(Order);
VALID_ORDER = logical(ismember(Order,[1,3,4,6,7,9,12,13,19,28,37]));

% IF desired order is not present,
if ~VALID_ORDER
    % THEN, set a default value
    Order = 7;
end

% BEGIN: define volume integration points and weights

switch Order
    
    case  1 % O(h), degree of precision 1
        Points  = [(1/3), (1/3)];
        Weights = (1/2);
    case  3 % O(h^3), degree of precision 2
        Points  = [(2/3), (1/6); (1/6), (2/3); (1/6), (1/6)];
        Weights = (1/2) * [1; 1; 1] / 3;
    case  4 % O(h^4), degree of precision 3
        Points  = [(1/3), (1/3); (1/5), (1/5); (3/5), (1/5); (1/5), (3/5)];
        Weights = (1/2) * [-27; 25; 25; 25] / 48;
    case  6 % O(h^6), degree of precision 4
        Points  = [0.816847572980459,  0.091576213509771;
                   0.091576213509771,  0.816847572980459;
                   0.091576213509771,  0.091576213509771;
                   0.108103018168070,  0.445948490915965;
                   0.445948490915965,  0.108103018168070;
                   0.445948490915965,  0.445948490915965];
        Weights = [0.109951743655322;
                   0.109951743655322;
                   0.109951743655322;
                   0.223381589678011;
                   0.223381589678011;
                   0.223381589678011] * (1/2);
    case  7 % O(h^7), degree of precision 5
        Points = [1./3.,                   1./3.;
                 (6.+sqrt(15.))    / 21., (6.+sqrt(15.))    / 21.;
                 (9.-2.*sqrt(15.)) / 21., (6.+sqrt(15.))    / 21.;
                 (6.+sqrt(15.))    / 21., (9.-2.*sqrt(15.)) / 21.;
                 (6.-sqrt(15.))    / 21., (6.-sqrt(15.))    / 21.;
                 (9.+2.*sqrt(15.)) / 21., (6.-sqrt(15.))    / 21.;
                 (6.-sqrt(15.))    / 21., (9.+2.*sqrt(15.)) / 21.];
        Weights = [ 9./40.;
                   (155.+sqrt(15.))/1200.;
                   (155.+sqrt(15.))/1200.;
                   (155.+sqrt(15.))/1200.;
                   (155.-sqrt(15.))/1200.;
                   (155.-sqrt(15.))/1200.;
                   (155.-sqrt(15.))/1200.] * (1/2);
    case  9 % O(h^9), degree of precision 6
        Points  = [0.124949503233232,  0.437525248383384;
                   0.437525248383384,  0.124949503233232;
                   0.437525248383384,  0.437525248383384;
                   0.797112651860071,  0.165409927389841;
                   0.797112651860071,  0.037477420750088;
                   0.165409927389841,  0.797112651860071;
                   0.165409927389841,  0.037477420750088;
                   0.037477420750088,  0.797112651860071;
                   0.037477420750088,  0.165409927389841];
        Weights = [0.205950504760887;
                   0.205950504760887;
                   0.205950504760887;
                   0.063691414286223;
                   0.063691414286223;
                   0.063691414286223;
                   0.063691414286223;
                   0.063691414286223;
                   0.063691414286223] * (1/2);
    case 12 % O(h^12), degree of precision 6
        Points  = [0.873821971016996,  0.063089014491502;
                   0.063089014491502,  0.873821971016996;
                   0.063089014491502,  0.063089014491502;
                   0.501426509658179,  0.249286745170910;
                   0.249286745170910,  0.501426509658179;
                   0.249286745170910,  0.249286745170910;
                   0.636502499121399,  0.310352451033785;
                   0.636502499121399,  0.053145049844816;
                   0.310352451033785,  0.636502499121399;
                   0.310352451033785,  0.053145049844816;
                   0.053145049844816,  0.636502499121399;
                   0.053145049844816,  0.310352451033785];
        Weights = [0.050844906370207;
                   0.050844906370207;
                   0.050844906370207;
                   0.116786275726379;
                   0.116786275726379;
                   0.116786275726379;
                   0.082851075618374;
                   0.082851075618374;
                   0.082851075618374;
                   0.082851075618374;
                   0.082851075618374;
                   0.082851075618374] * (1/2);
    case 13 % O(h^13), degree of precision 7
        Points  = [0.333333333333333,  0.333333333333333;
                   0.479308067841923,  0.260345966079038;
                   0.260345966079038,  0.479308067841923;
                   0.260345966079038,  0.260345966079038;
                   0.869739794195568,  0.065130102902216;
                   0.065130102902216,  0.869739794195568;
                   0.065130102902216,  0.065130102902216;
                   0.638444188569809,  0.312865496004875;
                   0.638444188569809,  0.048690315425316;
                   0.312865496004875,  0.638444188569809;
                   0.312865496004875,  0.048690315425316;
                   0.048690315425316,  0.638444188569809;
                   0.048690315425316,  0.312865496004875];
        Weights = [-0.149570044467670;
                   0.175615257433204;
                   0.175615257433204;
                   0.175615257433204;
                   0.053347235608839;
                   0.053347235608839;
                   0.053347235608839;
                   0.077113760890257;
                   0.077113760890257;
                   0.077113760890257;
                   0.077113760890257;
                   0.077113760890257;
                   0.077113760890257] * (1/2);
    case 19 % O(h^19), degree of precision 9
        Points  = [0.3333333333333333,  0.3333333333333333;
                   0.0206349616025259,  0.4896825191987370;
                   0.4896825191987370,  0.0206349616025259;
                   0.4896825191987370,  0.4896825191987370;
                   0.1258208170141290,  0.4370895914929355;
                   0.4370895914929355,  0.1258208170141290;
                   0.4370895914929355,  0.4370895914929355;
                   0.6235929287619356,  0.1882035356190322;
                   0.1882035356190322,  0.6235929287619356;
                   0.1882035356190322,  0.1882035356190322;
                   0.9105409732110941,  0.0447295133944530;
                   0.0447295133944530,  0.9105409732110941;
                   0.0447295133944530,  0.0447295133944530;
                   0.7411985987844980,  0.0368384120547363;
                   0.7411985987844980,  0.2219628891607657;
                   0.0368384120547363,  0.7411985987844980;
                   0.0368384120547363,  0.2219628891607657;
                   0.2219628891607657,  0.7411985987844980;
                   0.2219628891607657,  0.0368384120547363];
        Weights = [0.09713579628279610;
                   0.03133470022713983;
                   0.03133470022713983;
                   0.03133470022713983;
                   0.07782754100477543;
                   0.07782754100477543;
                   0.07782754100477543;
                   0.07964773892720910;
                   0.07964773892720910;
                   0.07964773892720910;
                   0.02557767565869810;
                   0.02557767565869810;
                   0.02557767565869810;
                   0.04328353937728940;
                   0.04328353937728940;
                   0.04328353937728940;
                   0.04328353937728940;
                   0.04328353937728940;
                   0.04328353937728940] * (1/2);
    case 28 % O(h^28), degree of precision 11
        Points  = [0.33333333333333333,  0.333333333333333333;
                   0.9480217181434233,   0.02598914092828833;
                   0.02598914092828833,  0.9480217181434233;
                   0.02598914092828833,  0.02598914092828833;
                   0.8114249947041546,   0.09428750264792270;
                   0.09428750264792270,  0.8114249947041546;
                   0.09428750264792270,  0.09428750264792270;
                   0.01072644996557060,  0.4946367750172147;
                   0.4946367750172147,   0.01072644996557060;
                   0.4946367750172147,   0.4946367750172147;
                   0.5853132347709715,   0.2073433826145142;
                   0.2073433826145142,   0.5853132347709715;
                   0.2073433826145142,   0.2073433826145142;
                   0.1221843885990187,   0.4389078057004907;
                   0.4389078057004907,   0.1221843885990187;
                   0.4389078057004907,   0.4389078057004907;
                   0.6779376548825902,   0.04484167758913055;
                   0.6779376548825902,   0.27722066752827925;
                   0.04484167758913055,  0.6779376548825902;
                   0.04484167758913055,  0.27722066752827925;
                   0.27722066752827925,  0.6779376548825902;
                   0.27722066752827925,  0.04484167758913055;
                   0.8588702812826364,   0.00000000000000000;
                   0.8588702812826364,   0.1411297187173636;
                   0.0000000000000000,   0.8588702812826364;
                   0.0000000000000000,   0.1411297187173636;
                   0.1411297187173636,   0.8588702812826364;
                   0.1411297187173636,   0.0000000000000000];
        Weights = [0.08797730116222190;
                   0.008744311553736190;
                   0.008744311553736190;
                   0.008744311553736190;
                   0.03808157199393533;
                   0.03808157199393533;
                   0.03808157199393533;
                   0.01885544805613125;
                   0.01885544805613125;
                   0.01885544805613125;
                   0.07215969754474100;
                   0.07215969754474100;
                   0.07215969754474100;
                   0.06932913870553720;
                   0.06932913870553720;
                   0.06932913870553720;
                   0.04105631542928860;
                   0.04105631542928860;
                   0.04105631542928860;
                   0.04105631542928860;
                   0.04105631542928860;
                   0.04105631542928860;
                   0.007362383783300573;
                   0.007362383783300573;
                   0.007362383783300573;
                   0.007362383783300573;
                   0.007362383783300573;
                   0.007362383783300573] * (1/2);
    case 37 % O(h^37), degree of precision 13
        Points  = [0.333333333333333333333333333333,  0.333333333333333333333333333333;
                   0.950275662924105565450352089520,  0.024862168537947217274823955239;
                   0.024862168537947217274823955239,  0.950275662924105565450352089520;
                   0.024862168537947217274823955239,  0.024862168537947217274823955239;
                   0.171614914923835347556304795551,  0.414192542538082326221847602214;
                   0.414192542538082326221847602214,  0.171614914923835347556304795551;
                   0.414192542538082326221847602214,  0.414192542538082326221847602214;
                   0.539412243677190440263092985511,  0.230293878161404779868453507244;
                   0.230293878161404779868453507244,  0.539412243677190440263092985511;
                   0.230293878161404779868453507244,  0.230293878161404779868453507244;
                   0.772160036676532561750285570113,  0.113919981661733719124857214943;
                   0.113919981661733719124857214943,  0.772160036676532561750285570113;
                   0.113919981661733719124857214943,  0.113919981661733719124857214943;
                   0.009085399949835353883572964740,  0.495457300025082323058213517632;
                   0.495457300025082323058213517632,  0.009085399949835353883572964740;
                   0.495457300025082323058213517632,  0.495457300025082323058213517632;
                   0.062277290305886993497083640527,  0.468861354847056503251458179727;
                   0.468861354847056503251458179727,  0.062277290305886993497083640527;
                   0.468861354847056503251458179727,  0.468861354847056503251458179727;
                   0.022076289653624405142446876931,  0.851306504174348550389457672223;
                   0.022076289653624405142446876931,  0.126617206172027096933163647918;
                   0.851306504174348550389457672223,  0.022076289653624405142446876931;
                   0.851306504174348550389457672223,  0.126617206172027096933163647918;
                   0.126617206172027096933163647918,  0.022076289653624405142446876931;
                   0.126617206172027096933163647918,  0.851306504174348550389457672223;
                   0.018620522802520968955913511549,  0.689441970728591295496647976487;
                   0.018620522802520968955913511549,  0.291937506468887771754472382212;
                   0.689441970728591295496647976487,  0.018620522802520968955913511549;
                   0.689441970728591295496647976487,  0.291937506468887771754472382212;
                   0.291937506468887771754472382212,  0.018620522802520968955913511549;
                   0.291937506468887771754472382212,  0.689441970728591295496647976487;
                   0.096506481292159228736516560903,  0.635867859433872768286976979827;
                   0.096506481292159228736516560903,  0.267625659273967961282458816185;
                   0.635867859433872768286976979827,  0.096506481292159228736516560903;
                   0.635867859433872768286976979827,  0.267625659273967961282458816185;
                   0.267625659273967961282458816185,  0.096506481292159228736516560903;
                   0.267625659273967961282458816185,  0.635867859433872768286976979827];
        Weights = [0.051739766065744133555179145422;
                   0.008007799555564801597804123460;
                   0.008007799555564801597804123460;
                   0.008007799555564801597804123460;
                   0.046868898981821644823226732071;
                   0.046868898981821644823226732071;
                   0.046868898981821644823226732071;
                   0.046590940183976487960361770070;
                   0.046590940183976487960361770070;
                   0.046590940183976487960361770070;
                   0.031016943313796381407646220131;
                   0.031016943313796381407646220131;
                   0.031016943313796381407646220131;
                   0.010791612736631273623178240136;
                   0.010791612736631273623178240136;
                   0.010791612736631273623178240136;
                   0.032195534242431618819414482205;
                   0.032195534242431618819414482205;
                   0.032195534242431618819414482205;
                   0.015445834210701583817692900053;
                   0.015445834210701583817692900053;
                   0.015445834210701583817692900053;
                   0.015445834210701583817692900053;
                   0.015445834210701583817692900053;
                   0.015445834210701583817692900053;
                   0.017822989923178661888748319485;
                   0.017822989923178661888748319485;
                   0.017822989923178661888748319485;
                   0.017822989923178661888748319485;
                   0.017822989923178661888748319485;
                   0.017822989923178661888748319485;
                   0.037038683681384627918546472190;
                   0.037038683681384627918546472190;
                   0.037038683681384627918546472190;
                   0.037038683681384627918546472190;
                   0.037038683681384627918546472190;
                   0.037038683681384627918546472190] * (1/2);
    otherwise
        error('Quadrature rule not chosen!');
end

% END: define volume integration points and weights

% END %

Contact us