0

I am having difficulty with the error "is not a list of numbers with dimensions..." when using FindRoot (and other numerical routines in Mathematica) to solve equations numerically when the argument of FindRoot is a function of the variable being chosen to solve the equation.

Here is simple code illustrating the behavior:

ClearAll["Global`*"]
rV = Table[r[i], {i, 2}];
pV = Table[p[i], {i, 2}];
aV = Table[a[i], {i, 2}];
aV0 = {1, 2}; 
f[aV_, rV_, pV_] := Table[aV[[i]] + rV[[i]] - pV[[i]], {i, 2}]
pVStar[aV_, rV_] := 
  pV /. FindRoot[f[aV, rV, pV] == 0, {{p[1], 1}, {p[2], 1}}]
g[aV_, rV_, pV_] := 2 aV + rV - pV
h[aV_, rV_] := g[aV, rV, pVStar[aV, rV]]
FindRoot[Table[h[aV0, rV][[i]] == 0, {i, 2}], {{r[1], 1}, {r[2], 1}}]

Here is the error:

During evaluation of In[1]:= FindRoot::nlnum: The function value {0. +r[1.],1. +r[2.]}
is not a list of numbers with dimensions {2} at {p[1],p[2]} = {1.,1.}. >>

I tried putting ?(VectorQ[#,NumericQ]&) after the rV_ in the definition of pVStar, as suggested in another thread about applying ?NumericQ to vectors:

pVStar[aV_, rV_?(VectorQ[#,NumericQ]&)] := 
 pV /. FindRoot[f[aV, rV, pV] == 0, {{p[1], 1}, {p[2], 1}}]

That gives the error

FindRoot::nveq: "The number of equations does not match the number of variables in
FindRoot[Table[h[aV0,rV][[i]]==0,{i,2}],{{r[1],1},{r[2],1}}]. 

Thoughts?

Dan
  • 51
  • 1
  • 4

1 Answers1

2

This is a way:

rV = Array[r, 2];    pV = Array[p, 2];    aV = Array[a, 2];
aV0 = {1, 2};

f[aV_, rV_, pV_] := aV + rV - pV
pVStar[aV_, rV_] := pV /. FindRoot[f[aV, rV, pV] == 0, Thread[{pV, 1}]]
g[aV_, rV_, pV_] := 2 aV + rV - pV
h[aV_?(VectorQ[#, NumericQ] &), rV_?(VectorQ[#, NumericQ] &)] := 
                                                        g[aV, rV, pVStar[aV, rV]]
h[aV_, rV_?(VectorQ[#, ! NumericQ@# &] &)] := {1, 1}

However:

FindRoot[h[aV0, rV], Thread[{rV, 1}]]

FindRoot::jsing: Encountered a singular Jacobian at the point {r[1],r[2]} = {1.,1.}. Try perturbing the initial point(s). >> h[aV_, rV_?(VectorQ[#, ! NumericQ@# &] &)] := {1, 1}

Because effectively your system has no solutions:

Plot3D[h[aV0, {x, y}], {x, -1, 1}, {y, -1, 1}]

Mathematica graphics

Dr. belisarius
  • 115,881
  • 13
  • 203
  • 453