<html>
  <head>

    <meta http-equiv="content-type" content="text/html; charset=utf-8">
  </head>
  <body bgcolor="#FFFFFF" text="#000000">
    Hi<br>
    <br>
    I've made a small script to play around with Bessel functions of the
    first kind... but this is very basic and I'm wondering if there's a
    smarter way.<br>
    <br>
    Below is the script I made (the function + a small test which plots
    the result).<br>
    <br>
    Best regards,<br>
    Claus<br>
    <br>
    <pre style="font-family:Monospaced;font-style:normal;font-size:12.0;"><span style="color:rgb(100,174,100);font-style:italic;">// bessel_test.sce</span>

<span style="color:rgb(176,24,19);">function</span> <span style="color:rgb(131,67,16);font-weight:bold;">z</span><span style="color:rgb(92,92,92);">=</span><span style="color:rgb(0,0,0);text-decoration:underline;">Jn</span><span style="color:rgb(74,85,219);">(</span><span style="color:rgb(131,67,16);font-weight:bold;">n</span><span style="color:rgb(0,0,0);">, </span><span style="color:rgb(131,67,16);font-weight:bold;">x</span><span style="color:rgb(74,85,219);">)</span>
    <span style="color:rgb(100,174,100);font-style:italic;">// The following power series approximates the nth Bessel</span>
    <span style="color:rgb(100,174,100);font-style:italic;">// function of the first kind for each input x</span>
    <span style="color:rgb(100,174,100);font-style:italic;">// In acoustics x = 2ka, defines the iput frequency k = omega / c and size</span>
    <span style="color:rgb(100,174,100);font-style:italic;">// of the piston radiator</span>
    <span style="color:rgb(0,0,0);">powerseries</span> <span style="color:rgb(92,92,92);">=</span> <span style="color:rgb(188,143,143);">0</span><span style="color:rgb(0,0,0);">;</span>
    <span style="color:rgb(160,32,240);">for</span> <span style="color:rgb(0,0,0);">m</span><span style="color:rgb(92,92,92);">=</span><span style="color:rgb(188,143,143);">0</span><span style="color:rgb(255,170,0);">:</span><span style="color:rgb(188,143,143);">19</span> <span style="color:rgb(100,174,100);font-style:italic;">// actually it should be infinity, but 10 approximates OK ...</span>
        <span style="color:rgb(0,0,0);">powerseries_m</span> <span style="color:rgb(92,92,92);">=</span> <span style="color:rgb(74,85,219);">(</span><span style="color:rgb(74,85,219);">(</span><span style="color:rgb(92,92,92);">-</span><span style="color:rgb(188,143,143);">1</span><span style="color:rgb(74,85,219);">)</span><span style="color:rgb(92,92,92);">^</span><span style="color:rgb(0,0,0);">m</span> <span style="color:rgb(92,92,92);">/</span> <span style="color:rgb(74,85,219);">(</span><span style="color:rgb(174,92,176);text-decoration:underline;">factorial</span><span style="color:rgb(74,85,219);">(</span><span style="color:rgb(0,0,0);">m</span><span style="color:rgb(74,85,219);">)</span> <span style="color:rgb(92,92,92);">*</span> <span style="color:rgb(174,92,176);text-decoration:underline;">factorial</span><span style="color:rgb(74,85,219);">(</span><span style="color:rgb(0,0,0);">m</span> <span style="color:rgb(92,92,92);">+</span> <span style="color:rgb(131,67,16);font-
weight:bold;">n</span><span style="color:rgb(74,85,219);">)</span><span style="color:rgb(74,85,219);">)</span><span style="color:rgb(74,85,219);">)</span> <span style="color:rgb(92,92,92);">*</span> <span style="color:rgb(74,85,219);">(</span><span style="color:rgb(131,67,16);font-weight:bold;">x</span><span style="color:rgb(92,92,92);">/</span><span style="color:rgb(188,143,143);">2</span><span style="color:rgb(74,85,219);">)</span><span style="color:rgb(92,92,92);">^</span><span style="color:rgb(74,85,219);">(</span><span style="color:rgb(188,143,143);">2</span><span style="color:rgb(92,92,92);">*</span><span style="color:rgb(0,0,0);">m</span><span style="color:rgb(74,85,219);">)</span><span style="color:rgb(0,0,0);">;</span>
        <span style="color:rgb(0,0,0);">powerseries</span> <span style="color:rgb(92,92,92);">=</span> <span style="color:rgb(0,0,0);">powerseries</span> <span style="color:rgb(92,92,92);">+</span> <span style="color:rgb(0,0,0);">powerseries_m</span><span style="color:rgb(0,0,0);">;</span> <span style="color:rgb(100,174,100);font-style:italic;">// sum the powerseries</span>
    <span style="color:rgb(160,32,240);">end</span>
    <span style="color:rgb(131,67,16);font-weight:bold;">z</span> <span style="color:rgb(92,92,92);">=</span> <span style="color:rgb(74,85,219);">(</span><span style="color:rgb(74,85,219);">(</span><span style="color:rgb(131,67,16);font-weight:bold;">x</span><span style="color:rgb(92,92,92);">/</span><span style="color:rgb(188,143,143);">2</span><span style="color:rgb(74,85,219);">)</span><span style="color:rgb(92,92,92);">^</span><span style="color:rgb(131,67,16);font-weight:bold;">n</span><span style="color:rgb(74,85,219);">)</span> <span style="color:rgb(92,92,92);">.*</span> <span style="color:rgb(0,0,0);">powerseries</span><span style="color:rgb(0,0,0);">;</span>
<span style="color:rgb(176,24,19);">endfunction</span>

<span style="color:rgb(0,0,0);">x_array</span> <span style="color:rgb(92,92,92);">=</span> <span style="color:rgb(188,143,143);">0</span><span style="color:rgb(255,170,0);">:</span><span style="color:rgb(188,143,143);">0.1</span><span style="color:rgb(255,170,0);">:</span><span style="color:rgb(188,143,143);">9.9</span><span style="color:rgb(0,0,0);">;</span> <span style="color:rgb(100,174,100);font-style:italic;">// define 100 points on the x-axis</span>
                     <span style="color:rgb(100,174,100);font-style:italic;">// with 2ka (x) from 0 to 10, and with m = 0-19,</span>
                     <span style="color:rgb(100,174,100);font-style:italic;">// this approximation is reasonably good for 2ka </span><span style="color:rgb(100,174,100);font-style:italic;"><</span><span style="color:rgb(100,174,100);font-style:italic;"> 10</span>

<span style="color:rgb(0,0,0);">z_0_array</span> <span style="color:rgb(92,92,92);">=</span> <span style="color:rgb(0,0,0);text-decoration:underline;">Jn</span><span style="color:rgb(74,85,219);">(</span><span style="color:rgb(188,143,143);">0</span><span style="color:rgb(0,0,0);">,</span><span style="color:rgb(0,0,0);">x_array</span><span style="color:rgb(74,85,219);">)</span><span style="color:rgb(0,0,0);">;</span> <span style="color:rgb(100,174,100);font-style:italic;">// Calculate J0</span>

<span style="color:rgb(0,0,0);">z_1_array</span> <span style="color:rgb(92,92,92);">=</span> <span style="color:rgb(0,0,0);text-decoration:underline;">Jn</span><span style="color:rgb(74,85,219);">(</span><span style="color:rgb(188,143,143);">1</span><span style="color:rgb(0,0,0);">,</span><span style="color:rgb(0,0,0);">x_array</span><span style="color:rgb(74,85,219);">)</span><span style="color:rgb(0,0,0);">;</span> <span style="color:rgb(100,174,100);font-style:italic;">// Calculate J1</span>

<span style="color:rgb(174,92,176);text-decoration:underline;">scf</span><span style="color:rgb(74,85,219);">(</span><span style="color:rgb(74,85,219);">)</span><span style="color:rgb(0,0,0);">;</span>
<span style="color:rgb(0,0,0);">a</span> <span style="color:rgb(92,92,92);">=</span> <span style="color:rgb(174,92,176);text-decoration:underline;">gca</span><span style="color:rgb(74,85,219);">(</span><span style="color:rgb(74,85,219);">)</span><span style="color:rgb(0,0,0);">;</span>
<span style="color:rgb(174,92,176);text-decoration:underline;">plot</span><span style="color:rgb(74,85,219);">(</span><span style="color:rgb(0,0,0);">x_array</span><span style="color:rgb(0,0,0);">,</span><span style="color:rgb(0,0,0);">z_0_array</span><span style="color:rgb(0,0,0);">,</span><span style="color:rgb(188,143,143);">'</span><span style="color:rgb(188,143,143);">-b</span><span style="color:rgb(188,143,143);">'</span><span style="color:rgb(74,85,219);">)</span><span style="color:rgb(0,0,0);">;</span>
<span style="color:rgb(174,92,176);text-decoration:underline;">plot</span><span style="color:rgb(74,85,219);">(</span><span style="color:rgb(0,0,0);">x_array</span><span style="color:rgb(0,0,0);">,</span><span style="color:rgb(0,0,0);">z_1_array</span><span style="color:rgb(0,0,0);">,</span><span style="color:rgb(188,143,143);">'</span><span style="color:rgb(188,143,143);">-r</span><span style="color:rgb(188,143,143);">'</span><span style="color:rgb(74,85,219);">)</span><span style="color:rgb(0,0,0);">;</span>
<span style="color:rgb(50,185,185);">xtitle</span><span style="color:rgb(74,85,219);">(</span><span style="color:rgb(188,143,143);">"</span><span style="color:rgb(188,143,143);">Bessel functions</span><span style="color:rgb(188,143,143);">"</span><span style="color:rgb(0,0,0);">,</span><span style="color:rgb(188,143,143);">"</span><span style="color:rgb(188,143,143);">x (2ka)</span><span style="color:rgb(188,143,143);">"</span><span style="color:rgb(0,0,0);">,</span><span style="color:rgb(188,143,143);">"</span><span style="color:rgb(188,143,143);">output (z)</span><span style="color:rgb(188,143,143);">"</span><span style="color:rgb(74,85,219);">)</span><span style="color:rgb(0,0,0);">;</span>
<span style="color:rgb(174,92,176);text-decoration:underline;">legend</span><span style="color:rgb(74,85,219);">(</span><span style="color:rgb(188,143,143);">"</span><span style="color:rgb(188,143,143);">J0</span><span style="color:rgb(188,143,143);">"</span><span style="color:rgb(0,0,0);">,</span><span style="color:rgb(188,143,143);">"</span><span style="color:rgb(188,143,143);">J1</span><span style="color:rgb(188,143,143);">"</span><span style="color:rgb(74,85,219);">)</span><span style="color:rgb(0,0,0);">;</span></pre>
    <br>
  </body>
</html>