無限区間数値積分

 投稿者：しばっち  投稿日：2019年 8月22日(木)20時32分43秒
  無限区間[-∞,∞]の数値積分を行います。

https://ja.wikipedia.org/wiki/数値積分

下記のように変数変換すると

  /∞           /1
  | f[x]dx ---> | {f[(1-t)/t]+f[-(1-t)/t]}/t^2dt
  /-∞          /0

積分区間は[0,1]となります。


被積分関数 f(x)=exp(-x*x)として

無限区間積分をリーマン和、台形則、シンプソン則等で行います。

  /∞
  | exp(-x*x)dx = sqr(π)
  /-∞

DECLARE EXTERNAL FUNCTION LEGENDRE,HERMITE,DE
PRINT INTEGRAL(8192)
PRINT INTEGRAL2(1024)
PRINT INTEGRAL3(1024)
PRINT LEGENDRE
PRINT HERMITE
PRINT DE
PRINT "----------------------"
PRINT SQR(PI) !'真値
END

EXTERNAL  FUNCTION FUNC(X)
LET FUNC=EXP(-X*X)
END FUNCTION

EXTERNAL FUNCTION INTEGRAL(N) !'リーマン和
LET A=0
LET B=1
LET H=(B-A)/(N+1)
FOR J=0 TO N
   LET T=A+H*J
   IF T<>0 THEN  LET S=S+H*(FUNC((1-T)/T)+FUNC(-(1-T)/T))/T/T
NEXT J
LET INTEGRAL=S
END FUNCTION

EXTERNAL  FUNCTION INTEGRAL2(N) !'台形則
DIM R(0 TO 1)
LET A=0
LET B=1
LET R(0)=1/2
LET R(1)=1/2
LET H=(B-A)/N
FOR K=0 TO N-1
   FOR J=0 TO 1
      LET T=A+H*(K+J)
      IF T<>0 THEN
         LET S=S+H*R(J)*(FUNC((1-T)/T)+FUNC(-(1-T)/T))/T/T
      END IF
   NEXT J
NEXT K
LET INTEGRAL2=S
END FUNCTION

EXTERNAL  FUNCTION INTEGRAL3(N) !'シンプソン則
DIM R(0 TO 2)
LET A=0
LET B=1
LET R(0)=1/3
LET R(1)=4/3
LET R(2)=1/3
LET H=(B-A)/N/2
FOR K=0 TO N-1
   FOR J=0 TO 2
      LET T=A+H*(2*K+J)
      IF T<>0 THEN
         LET S=S+H*R(J)*(FUNC((1-T)/T)+FUNC(-(1-T)/T))/T/T
      END IF
   NEXT J
NEXT K
LET INTEGRAL3=S
END FUNCTION

EXTERNAL  FUNCTION DE !'二重指数関数法(DE法)
LET H=1/1024
FOR T=-4 TO 4 STEP H
   LET S=S+H*FUNC(Q(T))*QQ(T)
NEXT T
LET DE=S

FUNCTION Q(X)
LET Q=SINH(PI/2*SINH(X))
END FUNCTION

FUNCTION QQ(X)
LET QQ=PI/2*COSH(X)*COSH(PI/2*SINH(X))
END FUNCTION
END FUNCTION

EXTERNAL  FUNCTION LEGENDRE !'ガウス・ルジャンドル則(100次)
LET A=0
LET B=1
LET U=(B+A)/2
LET V=(B-A)/2
FOR I=1 TO 100
   READ X,W
   LET T=U+V*X
   LET S=S+W*(FUNC((1-T)/T)+FUNC(-(1-T)/T))/T/T*V
NEXT I
LET LEGENDRE=S
DATA -.9997137267734412336782285,7.3463449050567173040632066E-04
DATA -.9984919506395958184001634,1.7093926535181052395293584E-03
DATA -.9962951347331251491861317,2.6839253715534824194395904E-03
DATA -.9931249370374434596520099,3.6559612013263751823424587E-03
DATA -.9889843952429917480044187,4.6244500634221193510957891E-03
DATA -.9838775407060570154961002,5.5884280038655151572119463E-03
DATA -.9778093584869182885537811,6.5469484508453227641521033E-03
DATA -.9707857757637063319308979,7.4990732554647115788287440E-03
DATA -.9628136542558155272936593,8.4438714696689714026208349E-03
DATA -.9539007829254917428493369,9.3804196536944579514182377E-03
DATA -.9440558701362559779627747,1.0307802574868969585782102E-02
DATA -.9332885350430795459243337,1.1225114023185977117221573E-02
DATA -.9216092981453339526669513,1.2131457662979497407744792E-02
DATA -.9090295709825296904671263,1.3025947892971542285558584E-02
DATA -.8955616449707269866985210,1.3907710703718772687954149E-02
DATA -.8812186793850184155733168,1.4775884527441301768879988E-02
DATA -.8660146884971646234107400,1.5629621077546002723936866E-02
DATA -.8499645278795912842933626,1.6468086176145212643104980E-02
DATA -.8330838798884008235429158,1.7290460568323582439344198E-02
DATA -.8153892383391762543939888,1.8095940722128116664390751E-02
DATA -.7968978923903144763895729,1.8883739613374904552941166E-02
DATA -.7776279096494954756275514,1.9653087494435305865381470E-02
DATA -.7575981185197071760356680,2.0403232646209432766838852E-02
DATA -.7368280898020207055124277,2.1133442112527641542672300E-02
DATA -.7153381175730564464599671,2.1843002416247386313953741E-02
DATA -.6931491993558019659486479,2.2531220256336272701796971E-02
DATA -.6702830156031410158025870,2.3197423185254121622488854E-02
DATA -.6467619085141292798326303,2.3840960265968205962560412E-02
DATA -.6226088602037077716041908,2.4461202707957052719975023E-02
DATA -.5978474702471787212648065,2.5057544481579589703764226E-02
DATA -.5725019326213811913168704,2.5629402910208116075642010E-02
DATA -.5465970120650941674679943,2.6176219239545676342308742E-02
DATA -.5201580198817630566468157,2.6697459183570962660384664E-02
DATA -.4932107892081909335693088,2.7192613446576880136491568E-02
DATA -.4657816497733580422492166,2.7661198220792388294204156E-02
DATA -.4378974021720315131089780,2.8102755659101173317648330E-02
DATA -.4095852916783015425288684,2.8516854322395097990936763E-02
DATA -.3808729816246299567633625,2.8903089601125203134876228E-02
DATA -.3517885263724217209723438,2.9261084110638276620119023E-02
DATA -.3223603439005291517224766,2.9590488059912642511754511E-02
DATA -.2926171880384719647375559,2.9890979593332830916836807E-02
DATA -.2625881203715034791689293,3.0162265105169144919068682E-02
DATA -.2323024818449739696495100,3.0404079526454820016507860E-02
DATA -.2017898640957359972360489,3.0616186583980448496459443E-02
DATA -.1710800805386032748875324,3.0798379031152590427713903E-02
DATA -.1402031372361139732075146,3.0950478850490988234063463E-02
DATA -.1091892035800611150034260,3.1072337427566516587810170E-02
DATA -.0780685828134366366948174,3.1163835696209906783818321E-02
DATA -.0468716824215916316149239,3.1224884254849357732376499E-02
DATA -.0156289844215430828722167,3.1255423453863356947642474E-02
DATA  .0156289844215430828722167,3.1255423453863356947642474E-02
DATA  .0468716824215916316149239,3.1224884254849357732376499E-02
DATA  .0780685828134366366948174,3.1163835696209906783818321E-02
DATA  .1091892035800611150034260,3.1072337427566516587810170E-02
DATA  .1402031372361139732075146,3.0950478850490988234063463E-02
DATA  .1710800805386032748875324,3.0798379031152590427713903E-02
DATA  .2017898640957359972360489,3.0616186583980448496459443E-02
DATA  .2323024818449739696495100,3.0404079526454820016507860E-02
DATA  .2625881203715034791689293,3.0162265105169144919068682E-02
DATA  .2926171880384719647375559,2.9890979593332830916836807E-02
DATA  .3223603439005291517224766,2.9590488059912642511754511E-02
DATA  .3517885263724217209723438,2.9261084110638276620119023E-02
DATA  .3808729816246299567633625,2.8903089601125203134876228E-02
DATA  .4095852916783015425288684,2.8516854322395097990936763E-02
DATA  .4378974021720315131089780,2.8102755659101173317648330E-02
DATA  .4657816497733580422492166,2.7661198220792388294204156E-02
DATA  .4932107892081909335693088,2.7192613446576880136491568E-02
DATA  .5201580198817630566468157,2.6697459183570962660384664E-02
DATA  .5465970120650941674679943,2.6176219239545676342308742E-02
DATA  .5725019326213811913168704,2.5629402910208116075642010E-02
DATA  .5978474702471787212648065,2.5057544481579589703764226E-02
DATA  .6226088602037077716041908,2.4461202707957052719975023E-02
DATA  .6467619085141292798326303,2.3840960265968205962560412E-02
DATA  .6702830156031410158025870,2.3197423185254121622488854E-02
DATA  .6931491993558019659486479,2.2531220256336272701796971E-02
DATA  .7153381175730564464599671,2.1843002416247386313953741E-02
DATA  .7368280898020207055124277,2.1133442112527641542672300E-02
DATA  .7575981185197071760356680,2.0403232646209432766838852E-02
DATA  .7776279096494954756275514,1.9653087494435305865381470E-02
DATA  .7968978923903144763895729,1.8883739613374904552941166E-02
DATA  .8153892383391762543939888,1.8095940722128116664390751E-02
DATA  .8330838798884008235429158,1.7290460568323582439344198E-02
DATA  .8499645278795912842933626,1.6468086176145212643104980E-02
DATA  .8660146884971646234107400,1.5629621077546002723936866E-02
DATA  .8812186793850184155733168,1.4775884527441301768879988E-02
DATA  .8955616449707269866985210,1.3907710703718772687954149E-02
DATA  .9090295709825296904671263,1.3025947892971542285558584E-02
DATA  .9216092981453339526669513,1.2131457662979497407744792E-02
DATA  .9332885350430795459243337,1.1225114023185977117221573E-02
DATA  .9440558701362559779627747,1.0307802574868969585782102E-02
DATA  .9539007829254917428493369,9.3804196536944579514182377E-03
DATA  .9628136542558155272936593,8.4438714696689714026208349E-03
DATA  .9707857757637063319308979,7.4990732554647115788287440E-03
DATA  .9778093584869182885537811,6.5469484508453227641521033E-03
DATA  .9838775407060570154961002,5.5884280038655151572119463E-03
DATA  .9889843952429917480044187,4.6244500634221193510957891E-03
DATA  .9931249370374434596520099,3.6559612013263751823424587E-03
DATA  .9962951347331251491861317,2.6839253715534824194395904E-03
DATA  .9984919506395958184001634,1.7093926535181052395293584E-03
DATA  .9997137267734412336782285,7.3463449050567173040632066E-04
END FUNCTION

EXTERNAL  FUNCTION HERMITE !'ガウス・エルミート則(50次)
FOR I=1 TO 50
   READ X,W
   LET S=S+FUNC(X)*W
NEXT I
LET HERMITE=S
DATA   -9.18240695812931737e+00 , 7.61348691118076750e-01
DATA   -8.52277103091780419e+00 , 5.88605297377289840e-01
DATA   -7.97562236820563655e+00 , 5.13304797851540890e-01
DATA   -7.48640942986419427e+00 , 4.68326211942547493e-01
DATA   -7.03432350977061065e+00 , 4.37553282300047560e-01
DATA   -6.60864797385535901e+00 , 4.14838821059022497e-01
DATA   -6.20295251927467162e+00 , 3.97244943775846315e-01
DATA   -5.81299467542040606e+00 , 3.83161392196479355e-01
DATA   -5.43578608722494814e+00 , 3.71619771249915029e-01
DATA   -5.06911758491723503e+00 , 3.61997265091617451e-01
DATA   -4.71129366616904279e+00 , 3.53872469979461721e-01
DATA   -4.36097316045457866e+00 , 3.46948769599378382e-01
DATA   -4.01706817285813439e+00 , 3.41010727258747851e-01
DATA   -3.67867706251526928e+00 , 3.35897876837685552e-01
DATA   -3.34503831393789109e+00 , 3.31488254247521860e-01
DATA   -3.01549776957452242e+00 , 3.27687661354606420e-01
DATA   -2.68948470226774507e+00 , 3.24422448708866827e-01
DATA   -2.36649390429866383e+00 , 3.21634537962055980e-01
DATA   -2.04607196868640921e+00 , 3.19277915993180742e-01
DATA   -1.72780654751589856e+00 , 3.17316124349604326e-01
DATA   -1.41131775489830006e+00 , 3.15720440230093699e-01
DATA   -1.09625112895768164e+00 , 3.14468550882471922e-01
DATA   -7.82271729554606886e-01 , 3.13543589989641611e-01
DATA   -4.69059056678236086e-01 , 3.12933448052389891e-01
DATA   -1.56302546889468676e-01 , 3.12630298030359123e-01
DATA   1.56302546889468676e-01 , 3.12630298030359123e-01
DATA   4.69059056678236086e-01 , 3.12933448052389891e-01
DATA   7.82271729554606886e-01 , 3.13543589989641611e-01
DATA   1.09625112895768164e+00 , 3.14468550882471922e-01
DATA   1.41131775489830006e+00 , 3.15720440230093699e-01
DATA   1.72780654751589856e+00 , 3.17316124349604326e-01
DATA   2.04607196868640921e+00 , 3.19277915993180742e-01
DATA   2.36649390429866383e+00 , 3.21634537962055980e-01
DATA   2.68948470226774507e+00 , 3.24422448708866827e-01
DATA   3.01549776957452242e+00 , 3.27687661354606420e-01
DATA   3.34503831393789109e+00 , 3.31488254247521860e-01
DATA   3.67867706251526928e+00 , 3.35897876837685552e-01
DATA   4.01706817285813439e+00 , 3.41010727258747851e-01
DATA   4.36097316045457866e+00 , 3.46948769599378382e-01
DATA   4.71129366616904279e+00 , 3.53872469979461721e-01
DATA   5.06911758491723503e+00 , 3.61997265091617451e-01
DATA   5.43578608722494814e+00 , 3.71619771249915029e-01
DATA   5.81299467542040606e+00 , 3.83161392196479355e-01
DATA   6.20295251927467162e+00 , 3.97244943775846315e-01
DATA   6.60864797385535901e+00 , 4.14838821059022497e-01
DATA   7.03432350977061065e+00 , 4.37553282300047560e-01
DATA   7.48640942986419427e+00 , 4.68326211942547493e-01
DATA   7.97562236820563655e+00 , 5.13304797851540890e-01
DATA   8.52277103091780419e+00 , 5.88605297377289840e-01
DATA   9.18240695812931737e+00 , 7.61348691118076750e-01
END FUNCTION
 		 

戻る