From what I understand, you want to compute the Jacobian matrix of a function $f$ at $(x_1,x_2,\ldots,x_n)$ and evaluate it at a given point, say $(2,2,\dots,2)$. This seems to be the function you are interested in: f[{X__}] := {X} RotateLeft[{X}, 1] - 1; In[2]:= f[{x1,x2,x3}] Out[2]= {-1+x1 x2,-1+x2 x3,-1+x1 x3} Since the dimension may change, we want to use `SetDelay` or `:=` so that the Jacobian function recomputes the matrix every time it is called. J[{X__}] := D[f[{X}], {{X}}] In[4]:= J[{x1,x2,x3}] Out[4]= {{x2,x1,0},{0,x3,x2},{x3,0,x1}} If you want the Jacobian matrix of $f$ at the point $(2,3,4)$, you can do the substitution In[5]:= J[{x1,x2,x3}] /. Thread[{x1,x2,x3}->{2,3,4}] Out[5]= {{3,2,0},{0,4,3},{4,0,2}} I would also recommend that you avoid using subscripts. Here is one way to generate a vector of variables such as `{x1,x2,x3,x4,x5}` programatically: In[6]:= ToExpression[Table["x"~~ToString[k],{k,1,5}]] Out[6]= {x1,x2,x3,x4,x5} variableVector[n_] := ToExpression[Table["x" ~~ ToString[k], {k, 1, n}]]; Now we can write a function which computes the Jacobian at a given point as follows: Jvalue[{Y__}] := With[ {X=variableVector[Length[{Y}]]}, J[X] /. Thread[X->{Y}] ] In[9]:= Jvalue[{2,3,4}] Out[9]= {{3,2,0},{0,4,3},{4,0,2}} **Update** If all you need is the sparse array of the Jacobian matrix of $f$ evaluated at a vector, you can do the following: Clear[x]; xVariables = Map[x, Range[10000]]; Evaluate[xVariables] = RandomInteger[{1, 20}, 10000]; sparse = SparseArray[{{10000, 10000} -> x[1], {10000, 1} -> x[10000], {i_,i_} -> x[i + 1], ({i_, j_} /; j - i == 1) -> x[i]}, {10000, 10000}] Here I used `x[i]` as dummy variable following J.M's suggestion, and `xVariables` is set to `{x[1],x[2],...,x[10000]}`. The next line sets `x[i]` to the `i`-th element of the random vector. The last line is where I use the fact that the Jacobian matrix of $f$ has a really nice pattern. You can check that `sparse` does have the correct values, say by evaluating `sparse[[1,2]]` and `sparse[[1,3]]`.