Tell me more ×
Mathematica Stack Exchange is a question and answer site for users of Mathematica. It's 100% free, no registration required.
up vote 2 down vote favorite

I am having trouble with the UnitStep function as in the title. My problem is very simple, but I am not able to get a numerical result.

I have 

f1[y] = 1/(E^((-1 + y)^2/2)*Sqrt[2*Pi])

g1[y] = (1.0028877725946312*^6*UnitStep[-7.963235463105154 - y])/
   E^((-1 + y)^2/2) + (0.12147136083763578*UnitStep[-7.963235463105154 + y])/
      E^((-1 + y)^2/2) + 
         1.001393070562657*
            (0.3484061634773921*Sqrt[E^(-(-1 + y)^2/2)] + 
                0.3484061634773921*Sqrt[E^(-(1 + y)^2/2)])^2*
            (-UnitStep[-7.963235463105154 + y] + UnitStep[7.963235463105154 + y])

and I want to solve the problem

$$N\left[\frac{1}{2}\int_{-\infty}^{\infty}\left(\sqrt{f1[y]}-\sqrt{g1[y]}\right)^2dy\right]$$

However, I did not get any result although I waited for a long time. I can plot $g1$ without any problem as well as $f1$, but I can not calculate the simple integral.

share|improve this question
there is a syntax error in your definition. Try f1[y_] = instead of f1[y] = – Thies Heidecke Feb 27 at 20:37
@ThiesHeidecke You mean f1[y_]:=, but as long as the NIntegrate uses y, that actually doesn't matter. I suspect the real problem is trying to use N[Integrate[...]] instead of NIntegrate. – Xerxes Feb 27 at 20:39
@Xerxes: in this case Set and SetDelayed are both fine. But you have a point with the N[Integrate[...]] construct. – Thies Heidecke Feb 27 at 20:41
@ThiesHeidecke even if I change to the other syntax I still have the same problem. – Seyhmus Güngören Feb 27 at 20:41
1  
Ah, got it. Well, you have your answer. When you wrap N outside of an integral, it first tries to evaluate it symbolically and only after it realises it can't (or until it succeeds), it tries numerical. Using NIntegrate it is done numerically from the start. The symbolic attempt is what takes long – Rojo Feb 27 at 20:46
show 5 more comments

1 Answer

up vote 6 down vote accepted

Using your definitions (using the placeholder pattern f1[y_] instead of the absolute pattern f1[y] is usually a good idea if you want to use it as a function that works with numerical values, too. Also using := (SetDelayed) instead of = (Set) inserts the left hand side value y into the definition every time you use it, which is closer to the behavior you would expect from a function):

f1[y_] := 1/(E^((-1 + y)^2/2)*Sqrt[2*Pi])

g1[y_] := (1.0028877725946312*^6*UnitStep[-7.963235463105154 - y])/
   E^((-1 + y)^2/2) + (0.12147136083763578*
     UnitStep[-7.963235463105154 + y])/E^((-1 + y)^2/2) + 
  1.001393070562657*(0.3484061634773921*Sqrt[E^(-(-1 + y)^2/2)] + 
      0.3484061634773921*
       Sqrt[E^(-(1 + y)^2/2)])^2*(-UnitStep[-7.963235463105154 + y] + 
     UnitStep[7.963235463105154 + y])

you could compute the integral via

(1/2) NIntegrate[(Sqrt[f1[y]]-Sqrt[g1[y]])^2, {y,-\[Infinity], \[Infinity]}]
(* 0.10271 *)

Using NIntegrate can save a lot of time, since Mathematica then knows that you are interested in a numerical solution from the start and doesn't waste time trying to find an analytical solution (Thanks to Xerxes and Rojo for pointing that out).

share|improve this answer
great! So my mistake was using the symbols. I guess in this case mathematica is searching first the algebraic solution and then trying to convert it to numerical value. Thanks alot. – Seyhmus Güngören Feb 27 at 20:49
Unless you Clear[y] prior to running this code, it is not safe to use Set instead of SetDelayed. Consider y=banana;f[y_]=y;f[0.1]. (* banana *) – Xerxes Feb 27 at 20:51
good point, i'll change the code for robustness. – Thies Heidecke Feb 27 at 20:52

Your Answer

 
discard

By posting your answer, you agree to the privacy policy and terms of service.

Not the answer you're looking for? Browse other questions tagged or ask your own question.