FindRoot with vector functions

Question

I'm trying to solve a system of non-linear equations with FindRoot, and I get the answer, but also a

*partd* (part spec longer than depth of object) 

error. This throws me off and makes me think it's poor programming. Is there a way to code this such that I don't get this message?

Clear["Global`*"]
F[x_] := {x[[2]]*x[[3]], -x[[1]]*x[[3]], -0.51 x[[1]]*x[[2]]}
dt = 0.1;
xi = {0, 1, 1};
xj = {x[1], x[2], x[3]};
FindRoot[xj == xi + dt*F[xj], {xj, xi}]

Part::partd: "Part specification xj[[2]] is longer than depth of object."
Part::partd: "Part specification xj[[3]] is longer than depth of object."
Part::partd: "Part specification xj[[1]] is longer than depth of object."

{{x[1], x[2], x[3]} -> {0.0985269, 0.990196, 0.995024}}

Answer 1

Define the helper functions and variables like this:

F[x_] := {x[[2]] x[[3]], -x[[1]] x[[3]], -0.51 x[[1]] x[[2]]}
dt = 0.1;
xi = {0, 1, 1};

Do not define xj since you need to use that as a single symbol in FindRoot.

Now this equation, xj == xi + dt F[xj] makes sense if we substitute a concrete vector value for xj. If we don't the equation will still evaluate with xj being a symbol (i.e. treated as a scalar), and cause errors for several obvious reasons:

In[5]:= xj == xi + dt F[xj]

During evaluation of In[5]:= Part::partd: Part specification xj[[2]] is longer than depth of object. >>

During evaluation of In[5]:= Part::partd: Part specification xj[[3]] is longer than depth of object. >>

During evaluation of In[5]:= Part::partd: Part specification xj[[1]] is longer than depth of object. >>

During evaluation of In[5]:= General::stop: Further output of Part::partd will be suppressed during this calculation. >>

Out[5]= xj == {0.1 xj[[2]] xj[[3]], 1 - 0.1 xj[[1]] xj[[3]], 1 - 0.051 xj[[1]] xj[[2]]}

The solution is to ask FindRoot not to evaluate the equation before it substitutes a value for xj:

In[6]:= FindRoot[xj == xi + dt F[xj], {{xj, xi}}, Evaluated -> False]

Out[6]= {xj -> {0.0985269, 0.990196, 0.995024}}

Yes, unfortunately using FindRoot with vector-variables is a bit less convenient then one would like.

---

I just realized that the Evaluated option of FindRoot is not documented. Here's a general and robust solution to these types of problems, without using the Evaluated option:

• Use FindRoot with a single function, like this: FindRoot[fun[x], {x, x0}]. It will find the value of x for which fun[x] is zero (where zero can be a vector/tensor of zeros).

• Make sure you give a starting value x0, and that x0 is a numerical vector (matrix, tensor, etc.) of the correct length (and dimension).

• Make sure fun[x] does not evaluate if x does not have a (vector-)value. If you evaluate fun[x] on its own, and x is a symbol, it should return as-is. To achieve this, define fun as

  Clear[fun]
  fun[x_?VectorQ] := ...
  

  or as

  Clear[fun]     
  fun[x_ ? (VectorQ[#, NumericQ]&)] := ...
  

  to require not only that x be a vector but also that each element of that vector be a numerical quantity (i.e. not a symbol)

Let's apply this to your problem:

In[13]:= fun[x_?(VectorQ[#, NumericQ] &)] := F[x] dt + xi - x

In[14]:= FindRoot[fun[x], {x, xi}]    
Out[14]= {x -> {0.0985269, 0.990196, 0.995024}}

Further reading (not the same problem but essentially analogous):

• http://support.wolfram.com/kb/3820