# How to remove points in time vector when I am getting a blowup Fisher information at those points please?

1 view (last 30 days)
Avan Al-Saffar on 18 Dec 2014
Edited: dpb on 18 Dec 2014
The code is :
function RunlogOscilnumericfishfixedn0omega2
omega=1;
N0=1;
k = 10;
A = 1;
p0 = 0.1;
t=(0.01:0.1:1000);
n=length(t);
[t,p] = ode45(@logOscilnumeric2,t,p0,[],omega,k,N0);
P = @(T) interp1(t,p,T);
v = (N0.^2.*P(t).^2.*(k.*sin(omega.*t) - P(t)).^4);
for i = 1 : n
if abs(v(i)) < 1e-15
f_v(i) = 0;
else
f_v(i) = ( (A.*( (k.^2 .*N0.*sin(omega.*t(i))) - (3.*N0.*P(t(i)).*k.*sin(omega.*t(i))) +
(k.^2.*omega.*cos(omega.*t(i))) + (2.*N0.*P(t(i)).^2) ).^2 )./v(i) ) ;
end
end
f = @(t) f_v;
I = integral( f, 0.01,1000,'ArrayValued',true)./1000
figure(1)
plot(t,I)
xlabel('Time')
ylabel('Fisher Information')
P = @(T) interp1(t,p,T);
g = @(t) ( (A.*( (k.^2 .*N0.*sin(omega.*t)) - (3.*N0.*P(t).*k.*sin(omega.*t)) +
(k.^2.*omega.*cos(omega.*t)) + (2.*N0.*P(t).^2) ).^2 )./ (N0.^2.*P(t).^2.*(k.*sin(omega.*t) - P
(t)).^4) );
figure(2)
plot(t,g(t))
1;
% function dpdt = logOscilnumeric2(t,p,omega,k,N0)
% dpdt = N0*sin(omega*t)*p - (N0 *p.^2/k);
% end