<div dir="ltr"><div><div>Dear Claus!<br><br></div>For me it seems that using standard functions is always better, than writing them from scratch.<br><br></div>So I recommend to use built-in Scilab functions <a href="http://help.scilab.org/docs/5.5.1/en_US/bessel.html" target="_blank">besseli, besselj, besselk, bessely, besselh</a>.<br><div><div><br>With best regards,<br>
<span><span>maintainer</span></span> of<a href="http://atoms.scilab.org/toolboxes/Mathieu" target="_blank"> Mathieu functions toolbox for Scilab</a>,<br>
IEEE member, Ph.D., <br>
Nikolay Strelkov.<br></div></div><div class="gmail_extra"><br clear="all">
<br><div class="gmail_quote">2015-01-24 19:53 GMT+03:00 Claus Futtrup <span dir="ltr"><<a href="mailto:cfuttrup@gmail.com" target="_blank">cfuttrup@gmail.com</a>></span>:<br><blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex">
<div 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>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>
</div>
<br>_______________________________________________<br>
users mailing list<br>
<a href="mailto:users@lists.scilab.org" target="_blank">users@lists.scilab.org</a><br>
<a href="http://lists.scilab.org/mailman/listinfo/users" target="_blank">http://lists.scilab.org/mailman/listinfo/users</a><br>
<br></blockquote></div><br></div></div>