function [ xtab, weight ] = lobatto_set ( norder ) %% LOBATTO_SET sets abscissas and weights for Lobatto quadrature. % % Discussion: % % The integration interval is [ -1, 1 ]. % % The weight function w(x) = 1. % % The integral to approximate: % % Integral ( -1 <= X <= 1 ) F(X) dX % % The quadrature rule: % % Sum ( 1 <= I <= NORDER ) WEIGHT(I) * F ( XTAB(I) ) % % The quadrature rule will integrate exactly all polynomials up to % X**(2*NORDER-3). % % The Lobatto rule is distinguished by the fact that both endpoints % (-1 and 1) are always abscissas of the rule. % % Licensing: % % This code is distributed under the GNU LGPL license. % % Modified: % % 12 October 2005 % % Author: % % John Burkardt % % Reference: % % Milton Abramowitz and Irene Stegun, % Handbook of Mathematical Functions, % National Bureau of Standards, 1964. % % Arthur Stroud and Don Secrest, % Gaussian Quadrature Formulas, % Prentice Hall, 1966. % % Daniel Zwillinger, editor, % Standard Mathematical Tables and Formulae, % 30th Edition, % CRC Press, 1996. % % Parameters: % % Input, integer NORDER, the order of the rule. % NORDER must be between 2 and 20. % % Output, real XTAB(NORDER), the abscissas for the rule. % % Output, real WEIGHT(NORDER), the weights of the rule. % The weights are positive, symmetric and should sum to 2. % if ( norder == 2 ) xtab(1) = - 1.0; xtab(2) = 1.0; weight(1) = 1.0; weight(2) = 1.0; elseif ( norder == 3 ) xtab(1) = - 1.0; xtab(2) = 0.0; xtab(3) = 1.0; weight(1) = 1.0 / 3.0; weight(2) = 4.0 / 3.0; weight(3) = 1.0 / 3.0; elseif ( norder == 4 ) xtab(1) = - 1.0; xtab(2) = - 0.447213595499957939281834733746; xtab(3) = 0.447213595499957939281834733746; xtab(4) = 1.0; weight(1) = 1.0 / 6.0; weight(2) = 5.0 / 6.0; weight(3) = 5.0 / 6.0; weight(4) = 1.0 / 6.0; elseif ( norder == 5 ) xtab(1) = - 1.0; xtab(2) = - 0.654653670707977143798292456247; xtab(3) = 0.0; xtab(4) = 0.654653670707977143798292456247; xtab(5) = 1.0; weight(1) = 9.0 / 90.0; weight(2) = 49.0 / 90.0; weight(3) = 64.0 / 90.0; weight(4) = 49.0 / 90.0; weight(5) = 9.0 / 90.0; elseif ( norder == 6 ) xtab(1) = - 1.0; xtab(2) = - 0.765055323929464692851002973959; xtab(3) = - 0.285231516480645096314150994041; xtab(4) = 0.285231516480645096314150994041; xtab(5) = 0.765055323929464692851002973959; xtab(6) = 1.0; weight(1) = 0.066666666666666666666666666667; weight(2) = 0.378474956297846980316612808212; weight(3) = 0.554858377035486353016720525121; weight(4) = 0.554858377035486353016720525121; weight(5) = 0.378474956297846980316612808212; weight(6) = 0.066666666666666666666666666667; elseif ( norder == 7 ) xtab(1) = - 1.0; xtab(2) = - 0.830223896278566929872032213967; xtab(3) = - 0.468848793470714213803771881909; xtab(4) = 0.0; xtab(5) = 0.468848793470714213803771881909; xtab(6) = 0.830223896278566929872032213967; xtab(7) = 1.0; weight(1) = 0.476190476190476190476190476190E-01; weight(2) = 0.276826047361565948010700406290; weight(3) = 0.431745381209862623417871022281; weight(4) = 0.487619047619047619047619047619; weight(5) = 0.431745381209862623417871022281; weight(6) = 0.276826047361565948010700406290; weight(7) = 0.476190476190476190476190476190E-01; elseif ( norder == 8 ) xtab(1) = - 1.0; xtab(2) = - 0.871740148509606615337445761221; xtab(3) = - 0.591700181433142302144510731398; xtab(4) = - 0.209299217902478868768657260345; xtab(5) = 0.209299217902478868768657260345; xtab(6) = 0.591700181433142302144510731398; xtab(7) = 0.871740148509606615337445761221; xtab(8) = 1.0; weight(1) = 0.357142857142857142857142857143E-01; weight(2) = 0.210704227143506039382991065776; weight(3) = 0.341122692483504364764240677108; weight(4) = 0.412458794658703881567052971402; weight(5) = 0.412458794658703881567052971402; weight(6) = 0.341122692483504364764240677108; weight(7) = 0.210704227143506039382991065776; weight(8) = 0.357142857142857142857142857143E-01; elseif ( norder == 9 ) xtab(1) = - 1.0; xtab(2) = - 0.899757995411460157312345244418; xtab(3) = - 0.677186279510737753445885427091; xtab(4) = - 0.363117463826178158710752068709; xtab(5) = 0.0; xtab(6) = 0.363117463826178158710752068709; xtab(7) = 0.677186279510737753445885427091; xtab(8) = 0.899757995411460157312345244418; xtab(9) = 1.0; weight(1) = 0.277777777777777777777777777778E-01; weight(2) = 0.165495361560805525046339720029; weight(3) = 0.274538712500161735280705618579; weight(4) = 0.346428510973046345115131532140; weight(5) = 0.371519274376417233560090702948; weight(6) = 0.346428510973046345115131532140; weight(7) = 0.274538712500161735280705618579; weight(8) = 0.165495361560805525046339720029; weight(9) = 0.277777777777777777777777777778E-01; elseif ( norder == 10 ) xtab(1) = - 1.0; xtab(2) = - 0.919533908166458813828932660822; xtab(3) = - 0.738773865105505075003106174860; xtab(4) = - 0.477924949810444495661175092731; xtab(5) = - 0.165278957666387024626219765958; xtab(6) = 0.165278957666387024626219765958; xtab(7) = 0.477924949810444495661175092731; xtab(8) = 0.738773865105505075003106174860; xtab(9) = 0.919533908166458813828932660822; xtab(10) = 1.0; weight(1) = 0.222222222222222222222222222222E-01; weight(2) = 0.133305990851070111126227170755; weight(3) = 0.224889342063126452119457821731; weight(4) = 0.292042683679683757875582257374; weight(5) = 0.327539761183897456656510527917; weight(6) = 0.327539761183897456656510527917; weight(7) = 0.292042683679683757875582257374; weight(8) = 0.224889342063126452119457821731; weight(9) = 0.133305990851070111126227170755; weight(10) = 0.222222222222222222222222222222E-01; elseif ( norder == 11 ) xtab(1) = - 1.0; xtab(2) = - 0.934001430408059134332274136099; xtab(3) = - 0.784483473663144418622417816108; xtab(4) = - 0.565235326996205006470963969478; xtab(5) = - 0.295758135586939391431911515559; xtab(6) = 0.0; xtab(7) = 0.295758135586939391431911515559; xtab(8) = 0.565235326996205006470963969478; xtab(9) = 0.784483473663144418622417816108; xtab(10) = 0.934001430408059134332274136099; xtab(11) = 1.0; weight(1) = 0.181818181818181818181818181818E-01; weight(2) = 0.109612273266994864461403449580; weight(3) = 0.187169881780305204108141521899; weight(4) = 0.248048104264028314040084866422; weight(5) = 0.286879124779008088679222403332; weight(6) = 0.300217595455690693785931881170; weight(7) = 0.286879124779008088679222403332; weight(8) = 0.248048104264028314040084866422; weight(9) = 0.187169881780305204108141521899; weight(10) = 0.109612273266994864461403449580; weight(11) = 0.181818181818181818181818181818E-01; elseif ( norder == 12 ) xtab(1) = - 1.0; xtab(2) = - 0.944899272222882223407580138303; xtab(3) = - 0.819279321644006678348641581717; xtab(4) = - 0.632876153031869677662404854444; xtab(5) = - 0.399530940965348932264349791567; xtab(6) = - 0.136552932854927554864061855740; xtab(7) = 0.136552932854927554864061855740; xtab(8) = 0.399530940965348932264349791567; xtab(9) = 0.632876153031869677662404854444; xtab(10) = 0.819279321644006678348641581717; xtab(11) = 0.944899272222882223407580138303; xtab(12) = 1.0; weight(1) = 0.151515151515151515151515151515E-01; weight(2) = 0.916845174131961306683425941341E-01; weight(3) = 0.157974705564370115164671062700; weight(4) = 0.212508417761021145358302077367; weight(5) = 0.251275603199201280293244412148; weight(6) = 0.271405240910696177000288338500; weight(7) = 0.271405240910696177000288338500; weight(8) = 0.251275603199201280293244412148; weight(9) = 0.212508417761021145358302077367; weight(10) = 0.157974705564370115164671062700; weight(11) = 0.916845174131961306683425941341E-01; weight(12) = 0.151515151515151515151515151515E-01; elseif ( norder == 13 ) xtab(1) = - 1.0; xtab(2) = - 0.953309846642163911896905464755; xtab(3) = - 0.846347564651872316865925607099; xtab(4) = - 0.686188469081757426072759039566; xtab(5) = - 0.482909821091336201746937233637; xtab(6) = - 0.249286930106239992568673700374; xtab(7) = 0.0; xtab(8) = 0.249286930106239992568673700374; xtab(9) = 0.482909821091336201746937233637; xtab(10) = 0.686188469081757426072759039566; xtab(11) = 0.846347564651872316865925607099; xtab(12) = 0.953309846642163911896905464755; xtab(13) = 1.0; weight(1) = 0.128205128205128205128205128205E-01; weight(2) = 0.778016867468189277935889883331E-01; weight(3) = 0.134981926689608349119914762589; weight(4) = 0.183646865203550092007494258747; weight(5) = 0.220767793566110086085534008379; weight(6) = 0.244015790306676356458578148360; weight(7) = 0.251930849333446736044138641541; weight(8) = 0.244015790306676356458578148360; weight(9) = 0.220767793566110086085534008379; weight(10) = 0.183646865203550092007494258747; weight(11) = 0.134981926689608349119914762589; weight(12) = 0.778016867468189277935889883331E-01; weight(13) = 0.128205128205128205128205128205E-01; elseif ( norder == 14 ) xtab(1) = - 1.0; xtab(2) = - 0.959935045267260901355100162015; xtab(3) = - 0.867801053830347251000220202908; xtab(4) = - 0.728868599091326140584672400521; xtab(5) = - 0.550639402928647055316622705859; xtab(6) = - 0.342724013342712845043903403642; xtab(7) = - 0.116331868883703867658776709736; xtab(8) = 0.116331868883703867658776709736; xtab(9) = 0.342724013342712845043903403642; xtab(10) = 0.550639402928647055316622705859; xtab(11) = 0.728868599091326140584672400521; xtab(12) = 0.867801053830347251000220202908; xtab(13) = 0.959935045267260901355100162015; xtab(14) = 1.0; weight(1) = 0.109890109890109890109890109890E-01; weight(2) = 0.668372844976812846340706607461E-01; weight(3) = 0.116586655898711651540996670655; weight(4) = 0.160021851762952142412820997988; weight(5) = 0.194826149373416118640331778376; weight(6) = 0.219126253009770754871162523954; weight(7) = 0.231612794468457058889628357293; weight(8) = 0.231612794468457058889628357293; weight(9) = 0.219126253009770754871162523954; weight(10) = 0.194826149373416118640331778376; weight(11) = 0.160021851762952142412820997988; weight(12) = 0.116586655898711651540996670655; weight(13) = 0.668372844976812846340706607461E-01; weight(14) = 0.109890109890109890109890109890E-01; elseif ( norder == 15 ) xtab(1) = - 1.0; xtab(2) = - 0.965245926503838572795851392070; xtab(3) = - 0.885082044222976298825401631482; xtab(4) = - 0.763519689951815200704118475976; xtab(5) = - 0.606253205469845711123529938637; xtab(6) = - 0.420638054713672480921896938739; xtab(7) = - 0.215353955363794238225679446273; xtab(8) = 0.0; xtab(9) = 0.215353955363794238225679446273; xtab(10) = 0.420638054713672480921896938739; xtab(11) = 0.606253205469845711123529938637; xtab(12) = 0.763519689951815200704118475976; xtab(13) = 0.885082044222976298825401631482; xtab(14) = 0.965245926503838572795851392070; xtab(15) = 1.0; weight(1) = 0.952380952380952380952380952381E-02; weight(2) = 0.580298930286012490968805840253E-01; weight(3) = 0.101660070325718067603666170789; weight(4) = 0.140511699802428109460446805644; weight(5) = 0.172789647253600949052077099408; weight(6) = 0.196987235964613356092500346507; weight(7) = 0.211973585926820920127430076977; weight(8) = 0.217048116348815649514950214251; weight(9) = 0.211973585926820920127430076977; weight(10) = 0.196987235964613356092500346507; weight(11) = 0.172789647253600949052077099408; weight(12) = 0.140511699802428109460446805644; weight(13) = 0.101660070325718067603666170789; weight(14) = 0.580298930286012490968805840253E-01; weight(15) = 0.952380952380952380952380952381E-02; elseif ( norder == 16 ) xtab(1) = - 1.0; xtab(2) = - 0.969568046270217932952242738367; xtab(3) = - 0.899200533093472092994628261520; xtab(4) = - 0.792008291861815063931088270963; xtab(5) = - 0.652388702882493089467883219641; xtab(6) = - 0.486059421887137611781890785847; xtab(7) = - 0.299830468900763208098353454722; xtab(8) = - 0.101326273521949447843033005046; xtab(9) = 0.101326273521949447843033005046; xtab(10) = 0.299830468900763208098353454722; xtab(11) = 0.486059421887137611781890785847; xtab(12) = 0.652388702882493089467883219641; xtab(13) = 0.792008291861815063931088270963; xtab(14) = 0.899200533093472092994628261520; xtab(15) = 0.969568046270217932952242738367; xtab(16) = 1.0; weight(1) = 0.833333333333333333333333333333E-02; weight(2) = 0.508503610059199054032449195655E-01; weight(3) = 0.893936973259308009910520801661E-01; weight(4) = 0.124255382132514098349536332657; weight(5) = 0.154026980807164280815644940485; weight(6) = 0.177491913391704125301075669528; weight(7) = 0.193690023825203584316913598854; weight(8) = 0.201958308178229871489199125411; weight(9) = 0.201958308178229871489199125411; weight(10) = 0.193690023825203584316913598854; weight(11) = 0.177491913391704125301075669528; weight(12) = 0.154026980807164280815644940485; weight(13) = 0.124255382132514098349536332657; weight(14) = 0.893936973259308009910520801661E-01; weight(15) = 0.508503610059199054032449195655E-01; weight(16) = 0.833333333333333333333333333333E-02; elseif ( norder == 17 ) xtab(1) = - 1.0; xtab(2) = - 0.973132176631418314156979501874; xtab(3) = - 0.910879995915573595623802506398; xtab(4) = - 0.815696251221770307106750553238; xtab(5) = - 0.691028980627684705394919357372; xtab(6) = - 0.541385399330101539123733407504; xtab(7) = - 0.372174433565477041907234680735; xtab(8) = - 0.189511973518317388304263014753; xtab(9) = 0.0; xtab(10) = 0.189511973518317388304263014753; xtab(11) = 0.372174433565477041907234680735; xtab(12) = 0.541385399330101539123733407504; xtab(13) = 0.691028980627684705394919357372; xtab(14) = 0.815696251221770307106750553238; xtab(15) = 0.910879995915573595623802506398; xtab(16) = 0.973132176631418314156979501874; xtab(17) = 1.0; weight(1) = 0.735294117647058823529411764706E-02; weight(2) = 0.449219405432542096474009546232E-01; weight(3) = 0.791982705036871191902644299528E-01; weight(4) = 0.110592909007028161375772705220; weight(5) = 0.137987746201926559056201574954; weight(6) = 0.160394661997621539516328365865; weight(7) = 0.177004253515657870436945745363; weight(8) = 0.187216339677619235892088482861; weight(9) = 0.190661874753469433299407247028; weight(10) = 0.187216339677619235892088482861; weight(11) = 0.177004253515657870436945745363; weight(12) = 0.160394661997621539516328365865; weight(13) = 0.137987746201926559056201574954; weight(14) = 0.110592909007028161375772705220; weight(15) = 0.791982705036871191902644299528E-01; weight(16) = 0.449219405432542096474009546232E-01; weight(17) = 0.735294117647058823529411764706E-02; elseif ( norder == 18 ) xtab(1) = - 1.0; xtab(2) = - 0.976105557412198542864518924342; xtab(3) = - 0.920649185347533873837854625431; xtab(4) = - 0.835593535218090213713646362328; xtab(5) = - 0.723679329283242681306210365302; xtab(6) = - 0.588504834318661761173535893194; xtab(7) = - 0.434415036912123975342287136741; xtab(8) = - 0.266362652878280984167665332026; xtab(9) = - 0.897490934846521110226450100886E-01; xtab(10) = 0.897490934846521110226450100886E-01; xtab(11) = 0.266362652878280984167665332026; xtab(12) = 0.434415036912123975342287136741; xtab(13) = 0.588504834318661761173535893194; xtab(14) = 0.723679329283242681306210365302; xtab(15) = 0.835593535218090213713646362328; xtab(16) = 0.920649185347533873837854625431; xtab(17) = 0.976105557412198542864518924342; xtab(18) = 1.0; weight(1) = 0.653594771241830065359477124183E-02; weight(2) = 0.399706288109140661375991764101E-01; weight(3) = 0.706371668856336649992229601678E-01; weight(4) = 0.990162717175028023944236053187E-01; weight(5) = 0.124210533132967100263396358897; weight(6) = 0.145411961573802267983003210494; weight(7) = 0.161939517237602489264326706700; weight(8) = 0.173262109489456226010614403827; weight(9) = 0.179015863439703082293818806944; weight(10) = 0.179015863439703082293818806944; weight(11) = 0.173262109489456226010614403827; weight(12) = 0.161939517237602489264326706700; weight(13) = 0.145411961573802267983003210494; weight(14) = 0.124210533132967100263396358897; weight(15) = 0.990162717175028023944236053187E-01; weight(16) = 0.706371668856336649992229601678E-01; weight(17) = 0.399706288109140661375991764101E-01; weight(18) = 0.653594771241830065359477124183E-02; elseif ( norder == 19 ) xtab(1) = - 1.0; xtab(2) = - 0.978611766222080095152634063110; xtab(3) = - 0.928901528152586243717940258797; xtab(4) = - 0.852460577796646093085955970041; xtab(5) = - 0.751494202552613014163637489634; xtab(6) = - 0.628908137265220497766832306229; xtab(7) = - 0.488229285680713502777909637625; xtab(8) = - 0.333504847824498610298500103845; xtab(9) = - 0.169186023409281571375154153445; xtab(10) = 0.0; xtab(11) = 0.169186023409281571375154153445; xtab(12) = 0.333504847824498610298500103845; xtab(13) = 0.488229285680713502777909637625; xtab(14) = 0.628908137265220497766832306229; xtab(15) = 0.751494202552613014163637489634; xtab(16) = 0.852460577796646093085955970041; xtab(17) = 0.928901528152586243717940258797; xtab(18) = 0.978611766222080095152634063110; xtab(19) = 1.0; weight(1) = 0.584795321637426900584795321637E-02; weight(2) = 0.357933651861764771154255690351E-01; weight(3) = 0.633818917626297368516956904183E-01; weight(4) = 0.891317570992070844480087905562E-01; weight(5) = 0.112315341477305044070910015464; weight(6) = 0.132267280448750776926046733910; weight(7) = 0.148413942595938885009680643668; weight(8) = 0.160290924044061241979910968184; weight(9) = 0.167556584527142867270137277740; weight(10) = 0.170001919284827234644672715617; weight(11) = 0.167556584527142867270137277740; weight(12) = 0.160290924044061241979910968184; weight(13) = 0.148413942595938885009680643668; weight(14) = 0.132267280448750776926046733910; weight(15) = 0.112315341477305044070910015464; weight(16) = 0.891317570992070844480087905562E-01; weight(17) = 0.633818917626297368516956904183E-01; weight(18) = 0.357933651861764771154255690351E-01; weight(19) = 0.584795321637426900584795321637E-02; elseif ( norder == 20 ) xtab(1) = - 1.0; xtab(2) = - 0.980743704893914171925446438584; xtab(3) = - 0.935934498812665435716181584931; xtab(4) = - 0.866877978089950141309847214616; xtab(5) = - 0.775368260952055870414317527595; xtab(6) = - 0.663776402290311289846403322971; xtab(7) = - 0.534992864031886261648135961829; xtab(8) = - 0.392353183713909299386474703816; xtab(9) = - 0.239551705922986495182401356927; xtab(10) = - 0.805459372388218379759445181596E-01; xtab(11) = 0.805459372388218379759445181596E-01; xtab(12) = 0.239551705922986495182401356927; xtab(13) = 0.392353183713909299386474703816; xtab(14) = 0.534992864031886261648135961829; xtab(15) = 0.663776402290311289846403322971; xtab(16) = 0.775368260952055870414317527595; xtab(17) = 0.866877978089950141309847214616; xtab(18) = 0.935934498812665435716181584931; xtab(19) = 0.980743704893914171925446438584; xtab(20) = 1.0; weight(1) = 0.526315789473684210526315789474E-02; weight(2) = 0.322371231884889414916050281173E-01; weight(3) = 0.571818021275668260047536271732E-01; weight(4) = 0.806317639961196031447768461137E-01; weight(5) = 0.101991499699450815683781205733; weight(6) = 0.120709227628674725099429705002; weight(7) = 0.136300482358724184489780792989; weight(8) = 0.148361554070916825814713013734; weight(9) = 0.156580102647475487158169896794; weight(10) = 0.160743286387845749007726726449; weight(11) = 0.160743286387845749007726726449; weight(12) = 0.156580102647475487158169896794; weight(13) = 0.148361554070916825814713013734; weight(14) = 0.136300482358724184489780792989; weight(15) = 0.120709227628674725099429705002; weight(16) = 0.101991499699450815683781205733; weight(17) = 0.806317639961196031447768461137E-01; weight(18) = 0.571818021275668260047536271732E-01; weight(19) = 0.322371231884889414916050281173E-01; weight(20) = 0.526315789473684210526315789474E-02; else fprintf ( 1, '\n' ); fprintf ( 1, 'LOBATTO_SET - Fatal error!\n' ); fprintf ( 1, ' Illegal value of NORDER = %d\n', norder ); fprintf ( 1, ' Legal values are between 2 and 20.\n' ); error ( 'LOBATTO_SET - Fatal error!' ); end