The error messages come from Integrate
/NIntegrate
calls on functions with symbolic (non-numerical) parameters a
, b
, c
. These calls come before FindRoot
substitutes numerical values for these parameters. (Taking N
of Integrate
here is effectively the same as calling NIntegrate
.)
In this case the best approach it seems to me is to do the integral first, since it can be done.
obj = Integrate[
x^{0, 2, 4} Exp[-(a*x^4 + b*x^2 + c)], {x, 0, Infinity},
Assumptions -> a > 0 && b > 0 && c > 0]
(*
{1/4 Sqrt[b/a] E^(b^2/(8 a) - c) BesselK[1/4, b^2/(8 a)],
(1/(16 Sqrt[2] Sqrt[a^3 b])) *
E^(b^2/(8 a) - c) π (-b^2 BesselI[-(1/4), b^2/(8 a)] +
(4 a + b^2) BesselI[1/4, b^2/(8 a)] +
b^2 (-BesselI[3/4, b^2/(8 a)] + BesselI[5/4, b^2/(8 a)])),
(Sqrt[b] E^(b^2/(8 a) - c) ((2 a + b^2) BesselK[1/4, b^2/(8 a)] -
b^2 BesselK[3/4, b^2/(8 a)]))/(32 a^(5/2))}
*)
aa = FindRoot[obj - {1, 2, 10}, {a, .01}, {b, .15}, {c, 1.33}]
(* {a -> 0.0108063, b -> 0.141937, c -> 0.671843} *)
Checks:
obj /. aa
(* {1., 2., 10.} *)
NIntegrate[x^{0, 2, 4} Exp[-(a*x^4 + b*x^2 + c)] /. aa, {x, 0, Infinity}]
(* {1., 2., 10.} *)
Further explanation:
Integrating vector expressions can be tricky. See NIntegrate over a list of functions and linked questions.
The integration is done on each component separately. NIntegrate
must see the components in the argument expression. The following, which is the normal way to use NumericQ
in numeric solvers like NIntegrate
, does not work because NIntegrate
decides the integrand is not a List
and gets confused when the values are not numbers:
i1[x_?NumericQ, a_?NumericQ, b_?NumericQ, c_?NumericQ] :=
(x^{0, 2, 4})*f[x, a, b, c]
FindRoot[NIntegrate[i1[x, a, b, c], {x, -Infinity, Infinity}] - {1, 2, 10},
{a, .01}, {b, .15}, {c, 1.33}]
NIntegrate::inumr: The integrand i1[x,a,b,c] has evaluated to non-numerical values.... >>
The proper way is to wrap the function up for FindRoot
instead:
obj2[a_?NumericQ, b_?NumericQ, c_?NumericQ] :=
NIntegrate[(x^{0, 2, 4})*f[x, a, b, c], {x, -Infinity, Infinity}]
But there's still one pitfall. If we substitute this for the integration there's a problem:
FindRoot[obj2[a, b, c] - {1, 2, 10}, {a, .01}, {b, .15}, {c, 1.33}]
FindRoot::nveq: The number of equations does not match the number of variables in FindRoot[obj2[a,b,c]-{1,2,10},{a,0.01},{b,0.15},{c,1.33}]. >>
The problem is that the argument to FindRoot
is evaluated symbolically before obj2
evaluates. You get a vector of obj2
calls:
obj2[a, b, c] - {1, 2, 10}
(* {-1 + obj2[a, b, c], -2 + obj2[a, b, c], -10 + obj2[a, b, c]} *)
When FindRoot
substitutes values for a
, b
, and c
, this evaluates and you get a matrix of values.
obj2[a, b, c] - {1, 2, 10} /. {a -> 0.01, b -> 0.15, c -> 1.33}
(*
{{0.0273318, 1.04919, 9.31436},
{-0.972668, 0.0491918, 8.31436},
{-8.97267, -7.95081, 0.314357}}
*)
The proper way to use obj2
is to use ==
instead of -
:
FindRoot[obj2[a, b, c] == {1, 2, 10}, {a, .01}, {b, .15}, {c, 1.33}]
(* {a -> 0.0108063, b -> 0.141937, c -> 1.36499} *)