1
$\begingroup$

The goal is to minimise a functional in this case I want to find the input function, defined by its (many) values on a grid for which the functional attains its minimum. Below is a piece of code where I i) make a grid, ii) evaluate some function on the grid (that will serve as the initial guess) iii) define some functions to interpolate and integrate iv) try to minimise the objective function (called total energy). All user defined functions seem to work fine, but when I evaluate FindMinimum NINtegrate throws error:

Integrate::inumr: ... has been evaluated to non-real values.

I read somewhere (on stack) that I should define an array to let Mathematica know I only want it to look for Real values on the grid points. I looked at this related question: Minimizing a function of many coordinates. (Note: I have limited the grid to 5 points, but I will in fact need 40000 to get accurate integrated values.)

genrhoi = Function[{rho1, h}, N[Range[rho1, 0.3, h]]]; 
grid = genrhoi[0., 1./16.] + 10^(-6)
densityongrid = (Exp[-3.*Abs[#1]] & ) /@ grid
setupintegrand = Function[{functionongrid}, functionongrid*grid^2.]; 
Integrator = Function[interpolatedintegrand, 4.*Pi*NIntegrate[interpolatedintegrand[x], 
      {x, grid[[1]], grid[[-1]]}, WorkingPrecision -> 10]]; 
Interpolator = Function[interpolatethis, Interpolation[
     Transpose[{grid, interpolatethis}], InterpolationOrder -> 4]]; 
nuclearpotential = Function[{position, charge, grid}, 
    (-1.*(charge/EuclideanDistance[position, #1]) & ) /@ grid]; 
externalpotential = nuclearpotential[0., 7., grid]; 
VRHO = Function[density, Integrator[Interpolator[setupintegrand[
       Evaluate[externalpotential*density]]]]]; 
TF = Function[density, Integrator[Interpolator[setupintegrand[
       Evaluate[density^(5./3.)]]]]]; 
totalenergy := TF[#1] + VRHO[#1] & ; 
vars = Array[a, Length[grid], 1]
FindMinimum[totalenergy[vars], Transpose[{vars, densityongrid}], MaxIterations -> 5, 
  PrecisionGoal -> 4, Gradient -> Evaluate[externalpotential + vars^(2./3.)], 
  Method -> "QuasiNewton", EvaluationMonitor -> Print[totalenergy[vars], vars, 
    ListLinePlot[Transpose[{grid, vars}]]]]
$\endgroup$
2
  • 2
    $\begingroup$ Alex, please use the help and this post to edit your question so it is legible and your code can be copy/pasted by others. Additionally, it would be a good idea for you to try and isolate a minimal code sample that exhibits the problematic behavior and post that alone. It will increase your chances of getting quality answers here. $\endgroup$
    – MarcoB
    Sep 9, 2015 at 20:40
  • $\begingroup$ To add to @MarcoB's comment: "I couldn't get it to work" is not very helpful for us. Can you be specific about what didn't work? Did it return wrong results? Did it throw errors? Did it return nothing at all? Did one of the command return unevaluated? Etc. $\endgroup$
    – march
    Sep 10, 2015 at 4:56

0

Browse other questions tagged or ask your own question.