Why the if statement does not work when I am trying to avoid the singularity in my code ??

function RunlogisticOscilnumericalfisherfixedn0omega

omega=1;

N0=1;

k = 10;

A = 1;

p0 = 0.1;

tspan=(0:0.1:4);

[t,p] = ode45(@logisticOscilnumerical,tspan,p0,[],omega,k,N0);

P = @(T) interp1(t,p,T);

f = @(t)  (  ( A.*( ( N0.* (sin(omega.*t)).^2 .*(1-(2.*P(t)./k))+(omega.*cos(omega.*t) ) ).^2 ) ./( (N0.^2).*((sin(omega.*t)).^4).*((P(t)-(P(t).^2./k)).^2 )) ) ) ;

v = ( (N0.^2).*((sin(omega.*t)).^4).*((P(t)-(P(t).^2./k)).^2 )); if v <= 1e-100

      f = 0

else

      f = I4

end

                                                                                I1 = integral( f, 0.01,1,'ArrayValued',true)./1

I2 = integral( f, 0.01,2,'ArrayValued',true)./2

I3 = integral( f, 0.01,3,'ArrayValued',true)./3

I4 = integral( f, 0.01,4,'ArrayValued',true)./4

I=[I1,I2,I3,I4] ;

T=[1,2,3,4] ;

figure(2)

plot(T,I)

g = @(t) (  ( A.*( ( N0.* (sin(omega.*t)).^2 .*(1-(2.*P(t)./k))+(omega.*cos(omega.*t) ) ).^2 ) ./( (N0.^2).*((sin(omega.*t)).^4).*((P(t)-(P(t).^2./k)).^2 )) ) ) ;

figure(3)

plot(tspan,g(tspan))

1;

% function dpdt = logisticOscilnumerical(t,p,omega,k,N0) 

% dpdt = N0*sin(omega*t)*p*(1-p/k);

% end

Could anyone help me to avoid singularity at x = 3.14159 using if statement please??

[Geoff Hayes]

Avan - it isn't clear to me why you have that if statement after the calculation of the integrals. Why set f to something if you are not going to use it?

One reason why your condition is failing is because v is a vector of dimension 41x1 rather than a scalar. With that in mind, what do you want your condition to be (and is it really necessary)?

[Avan Al-Saffar]

Thank you for your reply.

I am calling f as a function_handle in my integration. What I am doing is finding I1,..I4 but I am getting a warning error as follow : Warning: Minimum step size reached near x = 3.14159. There may be a singularity, or the tolerances may be too tight for this problem.

So I am trying by using if statement to solve this problem but I know that there is something wrong with the if statement, could you help me please?

[Torsten]

I already told you in another thread that your function f has singularities at 0, Pi/omega, 2*Pi/omega etc which are not removable (sin(omega*t) term in the denominator). So your integrals simply do not exist.

You will get I1=I2=I3=I4=infinity or -infinity.

Best wishes

Torsten.

[Avan Al-Saffar]

Ok, I know that my problem is with the term(sin(omega*t) but what if I tried to use if statement to avoid integration at those specific values, What do you think? Or do you have another idea?

Regards

Avan

[Torsten]

The values of a function at a finite number of points do not influence the integral.

Changing the function in a neighborhood of the singularities will give you a finite value for the integrals, but this finite value will depend on the size of the neighborhood you chose.

So I think it is better to reconsider how you came up with the function f and if this f is appropriate for your purpose.

Best wishes

Torsten.