[Scilab-users] Numderivative and integrate

[Hermes]

He,
here I put a solution for differentiation
  in your case we obtain The duration for Numederivative = 1.18248
and for this other: The duration for Numederivative = 0.00125101
the FourMatDiff function could be improved.
Is there any function in Scilab to create band Matrices?

//http://mailinglists.scilab.org/Scilab-users-Numderivative-and-integrate-td4037333.html
mode(0)
 
clear
 
//PATH = get_absolute_file_path("test_signal.sce");
 
// here we generate a sin signal just for testing
Gamma_ampl = 10; // in [g]
Frequency = 100; // in [Hz]
Omega = 2 * %pi * Frequency;
Number_of_periods = 10;
 
t = [0 : (1/(20 * Frequency)) : (Number_of_periods / Frequency)];
Acceleration = Gamma_ampl *sin(Omega * t);
Signal = [t ; Acceleration];
 
scf()
plot(t,Acceleration,'g');
xtitle('Ref curve - accelerogram','t','$\gamma \mbox{[m/s$^{2}$]}$');
xs2png(gcf(),'Accelerogram.png');
 
// signal first derivation : y' = (Omega * Gamma_ampl)*cos(Omega * t)
Acceleration_prime = (Omega * Gamma_ampl) * cos(Omega * t);
 
// calculation of the jerk -> numerical derivative
function y=Gamma_function(x)
    y = interpln(Signal,x)
endfunction
 
n = size(t,'*');
Numl_Gamma_derrivative = zeros(2,n);
Numl_Gamma_derrivative(1,:) = t;

Y=Gamma_function(Numl_Gamma_derrivative(1,:));
//
function M2=FourMatDiff2(n)
 for j=1:n
    for i=1:n
    //M(i,j) = 0; // main diagonal
    if pmodulo(i-j,n)==1 then
        M2(i,j)=-2/3
          elseif pmodulo(i-j,n)==2 then
        M2(i,j)=1/12
          elseif pmodulo(i-j,n)==n-1 then
        M2(i,j)=2/3
          elseif pmodulo(i-j,n)==n-2 then
        M2(i,j)=-1/12
          elseif i-j==0 then 
        M2(i,j)=0
        end
    end
end
//disp( full(M2) );   
endfunction;
N=n;
block=8;
sel=4;
h= t(3)-t(2);
X=0:h:N;
MD=FourMatDiff2(block);
tic();
//for i = 1 : n
//   Numl_Gamma_derrivative(2,i) = numderivative(Gamma_function,
Numl_Gamma_derrivative(1,i),1e-3,order=4);
//end
for i=1:N-block+1
 yd2($+1)=[MD/h*Y(i:i+block-1)'](sel);
end


Numederivative_duration = toc();
 
scf()
plot(Numl_Gamma_derrivative(1,4:197)',yd2,'b');
plot(t,Acceleration_prime,'ro');
xtitle('Jerk - numerical derivation vs derivative formula','t','$\mbox{J
[m/s}^{3}\mbox{]}$');
legend('Numerical','Reference',opt=1); 
xs2png(gcf(),'Jerk.png');

// information's
printf("The duration for Numederivative = %g\n",Numederivative_duration);

Gracias



--
Sent from: http://mailinglists.scilab.org/Scilab-users-Mailing-Lists-Archives-f2602246.html