I have a problem with a correct transfer of values of parameters when NIntegrate (and some subsequent functions) are coupled with minimizing procedures in Mathematica. Namely, let's define an operator:

V[r_,R_] := If[r < R, 0., 1.]
H[r_, R_] := m* D[#, r] + V[r, R]*# &;

which acts on a function (the function below is only an example but their nested form and 'ifs', which are important to me, would be possibly crucial)

X[R_, b_?NumberQ] := -b / Sin[b*R];
F1[r_, R_, b_?NumberQ] := If[r<=R, Sin[b*r] / r, Exp[-X[R, b]*r] / r]
F[r_, R_, b_?NumericQ] := F1[r, R, b]
m = 1.0;

I try to minimize this numerical integral (the best method seems to be 'FindMinValue' because initial value of minimizing parameter can be specified by hand). The two commands below do calculate it but non-numerical values are indicated in messages; how to avoid them?

(*1*) FindMinValue[NIntegrate[F[r,2.5,b]*H[r, 2.5]F[r,2.5,b]],{r,1,2.5}], {b,1.0}] (* R=2.5 is set here by hand *)
(*2*) With[{R=2.5},FindMinValue[NIntegrate[F[r,R,b]*H[r, R][F[r,R,b]],{r,1,R}],{b, 1.0}]]

But the following does not return any messages, surprisingly. What is the difference in transfer of 'R' into 'FindMinValue' here and in those previous cases? I cannot see a general rule.

(*3*) Do[{R = 2.5*i,abc = FindMinValue[NIntegrate[F[r,R,b]*H[r, R][F[r,R,b]],{r,1,R}],{b,1.0}],Print[abc]},{i,2}]

I would be grateful for help.

  • $\begingroup$ 1) Your F function is useless here, since F1 is already defined to evaluate only on numerical input; indeed if you swap F1 for F in your code you obtain the exact same results. 2) You had a missing opening bracket in (1): H[r, 2.5] F[r, 2.5, b]] should be H[r, 2.5][F[r, 2.5, b]], i.e. F` is the argument to H. $\endgroup$
    – MarcoB
    Feb 10, 2016 at 15:52
  • $\begingroup$ More importantly: all three approaches return exactly the same messages on my machine, and the same results. $\endgroup$
    – MarcoB
    Feb 10, 2016 at 15:54

1 Answer 1

  1. Your F wrapper function is doing nothing for you in your code, so I removed it and replaced it with direct calls to F1.
  2. NIntegrate in the argument to FindMinValue should not be evaluated unless it is passed explicitly numerical arguments, so it is best to wrap it in a function protected by NumericQ (functiontominimize below). Since all the other functions you defined are called from this one, they don't really need to include NumericQ checks.
  3. The minimization was running into numerical trouble. FindMinValue suggested to use more numerical precision, so I substituted your approximate values with exact ones, and increased the WorkingPrecision.

    Clear[V, H, X, F1]
    V[r_, R_] := If[r < R, 0., 1.]
    H[r_, R_] := m*D[#, r] + V[r, R]*# &;
    X[R_, b_] := -b/Sin[b*R];
    F1[r_, R_, b_] := If[r <= R, Sin[b*r]/r, Exp[-X[R, b]*r]/r]
    m = 1;
    functiontominimize[r_, R_, b_?NumericQ] := 
         NIntegrate[F1[r, R, b]*H[r, R][F1[r, R, b]], {r, 1, R}]
    With[{R = 25/10}, 
         FindMinValue[functiontominimize[r, R, b], b, WorkingPrecision -> 20]
    (* Out: -0.47610029093305050729 *)

The code above returns no message and executes faster as well.

  • $\begingroup$ Please excuse for some mistakes I made in original commands. Lack of '[' and that extra F appeared because these commands are taken from my original code stripped out from secondary parts. You solution, MarcoB, really solves the problem. Thank you very much. $\endgroup$
    – Marek
    Feb 10, 2016 at 19:21
  • $\begingroup$ Could you advice me some literature (a regular book on Mathematica) which describes the programming subtleties in details? I am a self-taught in Mathematica and my teaching now consists in solving my real problems by making a trial and trying to correct mistakes by guess. $\endgroup$
    – Marek
    Feb 10, 2016 at 19:34
  • $\begingroup$ @Marek take a look at these posts and the references therein: Where can I find examples of good Mathematica programming practice?; What Mathematica book to buy?; and the perhaps more advanced David Wagner's “Power programming with Mathematica”. $\endgroup$
    – MarcoB
    Feb 10, 2016 at 19:46
  • $\begingroup$ Ok, thank you very much. $\endgroup$
    – Marek
    Feb 10, 2016 at 19:49

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

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