MATLAB: Problem of integral combined with fzero

fzerointegral

Hi everyone,
I have encountered one problem when I want to use fzero to solve a integral equation, the detail is below
    rdata = normrnd(mu,sigma,1,n);
    xbar = mean(rdata);
    sd = std(rdata);
    aLehat = ( max( (xbar-T)*d/Du,(T-xbar)*d/Dl )./d_star )^2 + sd^2/d_star^2;
    delta = (xbar - T)./sd;
    b1 = @(y) sqrt(2./y).*delta;
    if delta >=0
        bL = @(y,c0) sqrt(n).*(Dl./d).*sqrt((2.*(aLehat./c0)./(n.*aLehat).*(  (d./Du.*delta).^2+1  ))./y-1);
        bU = @(y,c0) sqrt(n).*(Du./d).*sqrt((2.*(aLehat./c0)./(n.*aLehat).*(  (d./Du.*delta).^2+1  ))./y-1);
    else
        bL = @(y,c0) sqrt(n).*(Dl./d).*sqrt((2.*(aLehat./c0)./(n.*aLehat).*(  (d./Du.*delta).^2+1  ))./y-1);
        bU = @(y,c0) sqrt(n).*(Du./d).*sqrt((2.*(aLehat./c0)./(n.*aLehat).*(  (d./Du.*delta).^2+1  ))./y-1);
    end
    fun = @(y,c0) (    1./( exp( gammaln(alpha_no) ).*y.^(alpha_no+1) )   ).*exp(-1./y).* (   normcdf( b1(y)+bL(y,c0) ) - normcdf( b1(y)-bU(y,c0) )  );
    c0 = fzero(@(c0) (1-alpha) - integral(@(y) fun(y,c0), 0,(2.*(aLehat./ c0)./(n.*aLehat).*(  (d./Du.*delta).^2+1  ))),   0.01);
The problem is last line because when I put c0 in the denominator of 2.*(aLehat./ c0) , it could not work. However, if I put the c0 in numerator, it can work but the result is not what I want.
Hope everyone can give me some advice, thanks. Error message is below:
_Error using erfc Input must be real and full.
Error in normcdf>localnormcdf (line 124) p(todo) = 0.5 * erfc(-z ./ sqrt(2));
Error in normcdf (line 46) [varargout{1:max(1,nargout)}] = localnormcdf(uflag,x,varargin{:});
Error in @(y,c0)(1./(exp(gammaln(alpha_no)).*y.^(alpha_no+1))).*exp(-1./y).*(normcdf(b1(y)+bL(y,c0))-normcdf(b1(y)-bU(y,c0)))
Error in @(y)fun(y,c0)
Error in integralCalc/iterateScalarValued (line 314) fx = FUN(t);
Error in integralCalc/vadapt (line 132) [q,errbnd] = iterateScalarValued(u,tinterval,pathlen);
Error in integralCalc (line 75) [q,errbnd] = vadapt(@AtoBInvTransform,interval);
Error in integral (line 88) Q = integralCalc(fun,a,b,opstruct);
Error in @(c0)(1-alpha)-integral(@(y)fun(y,c0),0,(2.*(aLehat./c0)./(n.*aLehat).*((d./Du.*delta).^2+1)))
Error in fzero (line 363) a = x – dx; fa = FunFcn(a,varargin{:});

Best Answer

Instead of specifying a point (0.01) in the call to "fzero", you should specify a search interval which excludes the point c0=0.
Something like
x0 = [1e-6 1];
x = fzero(fun,x0)
Best wishes
Torsten.
Related Question