scilab/modules/differential_equations/macros/intsplin.sci

   1 // Scilab ( http://www.scilab.org/ ) - This file is part of Scilab
   2 // Copyright (C) INRIA
   3 // Copyright (C) 2017 - Samuel GOUGEON
   4 //
   5 // Copyright (C) 2012 - 2016 - Scilab Enterprises
   6 //
   7 // This file is hereby licensed under the terms of the GNU GPL v2.0,
   8 // pursuant to article 5.3.4 of the CeCILL v.2.1.
   9 // This file was originally licensed under the terms of the CeCILL v2.1,
  10 // and continues to be available under such terms.
  11 // For more information, see the COPYING file which you should have received
  12 // along with this program.
  13
  14 function v = intsplin(x, s)
  15     //splin  numerical integration.
  16     //v = intsplin(x,s) computes the integral of y with respect to x using
  17     //splin interpolation and integration.
  18     //x and y must be vectors of the same dimension
  19     //
  20     //v = intsplin(s) computes the integral of y assuming unit
  21     //spacing between the data points.
  22
  23     x = x(:);
  24     if argn(2) < 2 then
  25         s = x;
  26         x = (1:size(s,"*"))';
  27     else
  28         s = s(:);
  29         if size(x,"*")<>size(s,"*") then
  30             msg = _("%s: Wrong size for input arguments: Same size expected.\n")
  31             error(msprintf(msg, "intsplin"));
  32         end
  33     end
  34     if ~isreal(x) then
  35         if isreal(x,%eps)
  36             x = real(x);
  37         else
  38             msg = _("%s: Argument #%d: Real number(s) expected.\n");
  39             error(msprintf(msg, "intsplin",1))
  40         end
  41     end
  42     h = x(2:$) - x(1:$-1);
  43     y = real(s);
  44     d = splin(x, y);
  45     v = sum((h.*(d(1:$-1)-d(2:$))/12 + (y(1:$-1)+y(2:$))/2).*h);
  46     if ~isreal(s) then
  47         y = imag(s);
  48         d = splin(x, y);
  49         v = v + %i*sum((h.*(d(1:$-1)-d(2:$))/12 + (y(1:$-1)+y(2:$))/2).*h);
  50     end
  51 endfunction