Mathematica Stack Exchange is a question and answer site for users of Mathematica. It's 100% free, no registration required.

Sign up
Here's how it works:
  1. Anybody can ask a question
  2. Anybody can answer
  3. The best answers are voted up and rise to the top

This post has been edited heavily in response to the comments: my initial question was very confusing.

I am trying to write my game-theoretic model in Wolfram Language. There is only one remaining step and I am not sure if it is within reach of Wolfram Language (or Mathematica). I describe a simplified version of my problem below.

Cournot game

There are $M$ firms. Each firm $i$ can produce some output $q_i\geq0$ at a cost of $c_i$ per unit. Once the firms have (simultaneously) decided on their outputs, the price is set in the market according to an inverse demand function $P(Q)=\alpha-Q$, where $Q=q_1+\ldots+q_M$ is the aggregate output and $\alpha>0$ is market size parameter. Firm $i$'s profits are $q_iP(Q)-q_ic_i=q_i(\alpha-q_1-\ldots-q_M-c_i)$. Each firm maximizes its profits taking outputs of the opponents as given. By differentiating the profit function with respect to $q_i$ and imposing the necessary condition for maximizing a function, I obtain the equilibrium outputs: $$q_i^*=\frac{1}{2}\left(\alpha-\sum_{k\neq i}q_k^*-c_k\right) \text{ for } i=\{1,\ldots,M\}.$$ A quick derivation on the piece of paper gives me the explicit solution: $$q_i^*=\frac{1}{M+1}\left(\alpha-Mc_i+\sum_{k\neq i}c_k\right)$$ Here's the thing: I would like to obtain exactly this result (the above equation) in Mathematica. Is it possible? Note: I don't want to substitute any numbers for $M, c_1, \ldots, c_M$ or $\alpha$. They are general, exogenous parameters.

In other words, the above problem boils down to inverting a $M\times M$ matrix of the form: $$\left[\begin{array}{cccc} 2 & 1 & \ldots & 1\\ 1 & 2 & \ldots & 1\\ \vdots & \vdots & \ddots & \vdots \\ 1 & 1 & \ldots & 2 \end{array}\right].$$ Can one do this without specifying $M$?

I hope I expressed myself clearer this time.

share|improve this question
If you have vectors c and q both of length m (any m) then Total[c q] does your calculation. You can create such vectors in many ways: ConstantArray[0,m] creates an m-length vector of zeros. RandomReal[{0,1},m] creates a vector of length m each of which is a random number in {0,1}. – bill s Apr 22 '14 at 17:48
No, that is not what I aim for. I have updated the description. I do not want to use specific values. – Krzysztofik Apr 22 '14 at 18:00
I feel the closest you could get to a vector of symbolic length would be SymbolicTensors. However, they are mostly used in vector calculus like stuff (proving identities and so) and may not be what you are looking for. Perhaps it would help if you could elaborate a bit on what you want to do with them. – Sjoerd C. de Vries Apr 22 '14 at 18:11
I extended the description, I hope that clarifies things a bit. – Krzysztofik Apr 22 '14 at 19:14
IdentityMatrix[M] - 1/(M + 1) ConstantArray[1, M]? – kglr Apr 23 '14 at 17:50

Maybe something like this:

cournotEq[m_] := (IdentityMatrix[m] - 
    1/(m + 1) ConstantArray[1, m]).Array[\[Alpha] - Subscript[c, #] &,  m]

$\left\{\frac{1}{3} \left(\alpha -2 c_1+c_2\right),\frac{1}{3} \left(\alpha +c_1-2 c_2\right)\right\}$


$\left\{\frac{1}{4} \left(\alpha -3 c_1+c_2+c_3\right),\frac{1}{4} \left(\alpha +c_1-3 c_2+c_3\right),\frac{1}{4} \left(\alpha +c_1+c_2-3 c_3\right)\right\} $

share|improve this answer
Almost. Can one make it work for cournotEq[M]? Without giving a concrete number? As it is, it is not possible. I want to get the exact equation I have given in the description (e.g. in a matrix form). – Krzysztofik Apr 23 '14 at 20:59
@KrzysztofBrzezinski, I don't think it is possible to have it work for symbolic M. – kglr Apr 23 '14 at 21:23

You can invert the matrix $$\left[\begin{array}{cccc} 2 & 1 & \ldots & 1\\ 1 & 2 & \ldots & 1\\ \vdots & \vdots & \ddots & \vdots \\ 1 & 1 & \ldots & 2 \end{array}\right].$$ (for general m) using the Woodbury matrix identity. First, let one be the m-vector of all ones. Then the matrix you want to invert is expressible as

mat = IdentityMatrix[m] + Outer[Times, one, one];

To use the matrix inversion lemma (with the notation from the wikipedia article), let A = IdentityMatrix[m], U=one, V=Transpose[one], C=A; and the lemma says the inverse is:

matInv= IdentityMatrix[m] - Outer[Times, one, one]/(m + 1);

This holds for all m. We can check this for an specific m like so:

m = 5;
one = ConstantArray[1, m];
eye = IdentityMatrix[m];
mat = eye + Outer[Times, one, one];
matInv = eye - Outer[Times, one, one]/(m + 1);
mat.matInv == matInv.mat == eye

to which the response is True. To do this truly within Mathematica, one approach would be to incorporate knowledge of the inversion lemma.

share|improve this answer
I'm afraid this does not solve the problem. Can one make this method work without specifying m? – Krzysztofik Apr 24 '14 at 17:10
If you want a numerical inverse to the matrix, then you must specify m. If you don't want a numerical inverse, then the formula for matInv is the answer: $I-1.1'/(m+1)$. – bill s Apr 24 '14 at 19:13

Your Answer


By posting your answer, you agree to the privacy policy and terms of service.

Not the answer you're looking for? Browse other questions tagged or ask your own question.