0

I am using matlab and enter three functions: Integral,AddtoIntegralTerm1 and AddtoIntegralTerm2. I want to add the three numbers, so I write it as follow Integral(0,139)+AddtoIntegralTerm1(0,139)+AddtoIntegralTerm2(0,139) which gives me a number 1048576.

However, when I divided it by 1e21 (I calculate the order of each number which is 1e21) and add, i.e. Integral(0,139)/1e21+AddtoIntegralTerm1(0,139)/1e21+AddtoIntegralTerm2(0,139)/1e21, I get a number 0. Is it related to precision? How to fix it? The functions are listed below:

        % Use of indefinite integral
    pit=@(t)-MinOrMax.*...
        (-1./lambda.*(exp(-lambda.*min(t,CutOff)))-(lambda.*beta)./(sigma.^2).*exp(-lambda.*t));
    pitau=@(tau)MinOrMax.*...
        (-1./lambda.*(exp(-lambda.*min(tau,CutOff)))-(lambda.*beta)./(sigma.^2).*exp(-lambda.*tau));
    vt=@(t)-MinOrMax.*sigma.^2.*...
        (1./(lambda.^2).*(exp(lambda.*min(t,CutOff)))+...
        ...
        exp(-lambda.*CutOff)./(2.*lambda.^2).*(exp(2.*lambda.*t)-exp(2.*lambda.*min(t,CutOff)))+...
        ...
        (lambda.*beta)./(sigma.^2.*lambda).*(exp(lambda.*t)));
    vtau=@(tau)-MinOrMax.*sigma.^2.*...
        (1./(lambda.^2).*(exp(lambda.*min(tau,CutOff)))+...
        ...
        exp(-lambda.*CutOff)./(2.*lambda.^2).*(exp(2.*lambda.*tau)-exp(2.*lambda.*min(tau,CutOff)))+...
        ...
        (lambda.*beta)./(sigma.^2.*lambda.^2).*(exp(lambda.*tau)));

    Integral1=...
        @(t,tau)sigma.^2./(2.*lambda.^2).*...
        (1./(2.*lambda).*exp(-2.*lambda.*CutOff).*(exp(2.*lambda.*tau)-exp(2.*lambda.*t)).*(CutOff<t)+...
        ...
        ((CutOff-t)+exp(-2.*lambda.*CutOff)./(2.*lambda).*...
        (exp(2.*lambda.*tau)-exp(2.*lambda.*CutOff))).*((CutOff<=tau).*(CutOff>=t))+...
        ...
        (tau-t).*(CutOff>=tau));
    Integral2=...
        @(t,tau)2.*((lambda.*beta)./sigma.^2).*sigma.^2./(2.*lambda).*...
        (1./lambda.*exp(-lambda.*CutOff).*(exp(lambda.*tau)-exp(lambda.*t)).*(CutOff<t)+...
        ...
        ((CutOff-t)+exp(-lambda.*CutOff)./(lambda).*...
        (exp(lambda.*tau)-exp(lambda.*CutOff))).*((CutOff<=tau).*(CutOff>=t))+...
        ...
        (tau-t).*(CutOff>=tau));
    Integral3=@(t,tau)sigma.^2./2.*((lambda.*beta)./sigma.^2).^2.*(tau-t);
    Integral=@(t,tau)Integral1(t,tau)+Integral2(t,tau)+Integral3(t,tau);

    AddtoIntegralTerm1=@(t,tau)sigma.^2./(2.*lambda).*pitau(tau).^2.*(exp(2.*lambda.*tau)-exp(2.*lambda.*t))./2;
    AddtoIntegralTerm2=@(t,tau)-pitau(tau).*vtau(tau);

Thanks for the help.

Ander Biguri
  • 35,140
  • 11
  • 74
  • 120
will_cheuk
  • 379
  • 3
  • 12
  • Yes its related to precision. Your only options is `vpa`. Otherwise, this problem is not related to MATLAB, but to computers in general. – Ander Biguri Aug 08 '18 at 12:35
  • Thanks. What should I do to fix it then? I have just checked that vpa is for symbolic function but I get anonymous function. What should I code to solve it? – will_cheuk Aug 08 '18 at 12:43
  • The first thing you should definitely do is think why would you want to divide a number by such a massive other number, and try to think alternatives. Why do you want to do that? – Ander Biguri Aug 08 '18 at 12:45
  • Originally I do not want to divide a number. I want to use the value calculated to evaluate exponential function. For example, I want to compute when t=0, tau=139, CutOff=0.05, what is the exp(Integral(0,139)+AddtoIntegralTerm1(0,139)+AddtoIntegralTerm2(0,139)). However, when I compute manually after some simplification, I find that all the terms which make Integral(0,139)+AddtoIntegralTerm1(0,139)+AddtoIntegralTerm2(0,139) to be very large can be cancelled out to give a small number so that its exponential can be computed. However, if I directly sub this, I get a very large number and test. – will_cheuk Aug 08 '18 at 12:54
  • You sayd that the sum of those things is `1048576`. That is not "cancelled". – Ander Biguri Aug 08 '18 at 12:56
  • Yes, when I use matlab directly without division, I get 1048576. However, when I do simplification to the formula after substituting t=0, tau=139, CutOff=0.05 manually first, those terms which give very large numbers disappear. So I suspect it relates to precision of number and check for it. Since each number I compute is 1e21, so I divide 1e21 and add it, it gives 0. I also do the division from 1e1 to 1e20 as well, in most cases, it gives 0 (some are not, say 1e3). I wonder if it is related to precision and ask how to fix it since I would like to use the same code for other uses myself. – will_cheuk Aug 08 '18 at 13:07
  • Yes, its related to precision. Read here: https://stackoverflow.com/questions/686439/why-is-24-0000-not-equal-to-24-0000-in-matlab – Ander Biguri Aug 08 '18 at 13:12
  • Thanks. Since later, I would like to pass out a function g=@(t,tau)exp(Integral(t,tau)+AddtoIntegralTerm1(t,tau)+AddtoIntegralTerm2(t,tau)) to another script and do integration (integral(@(tau)g(0,tau),0,200)), the answer gives inf/NaN, while it should be finite, how to fix it? It is impossible for me to further simplify it myself, my simplification is correct when tau>CutOff, otherwise, we still need to use the same function displayed. – will_cheuk Aug 08 '18 at 13:21

0 Answers0