**Competing in different settings**

In this project we consider 2 firms who compete in the same duopolistic market. We will look at three possible competition forms, which are characterized by

**Cournot**

- Firms compete in quantities, and decide upon these independently and simultaneously 

- Firms profit maximize given the other's choice of quantity

- Both firms have market power, thus one's decision of output affects the price of the other's goods

**Bertrand**

- Firms compete in prices, and decide upon these simultaneously and as if they were in a perfect competition setting

- We assume that consumers seek the lowest price (assuming the two products are homogenous)

- When prices are set, output is given as 

$$D\left(p_{i},p_{j}\right)=\begin{cases}
D\left(p_{i}\right) & p_{i}<p_{j}\\
\frac{D\left(p_{i}\right)}{2} & p_{i}=p_{j}\\
0 & p_{i}>p_{j}
\end{cases}$$

**Stackelberg**

- Firms compete in quantities, and decide upon these sequentially. 

- One firm is assumed to be the leader, and the other is the follower. 

In the following we will analyse a sitaution with linear demand with $a=101$, $b=0.1$ and cost equal to 1 for all firms under all three situations.


```python
#Import packages
from scipy import optimize,arange
from numpy import array
import numpy as np
import matplotlib.pyplot as plt
import random
import sympy as sm
from math import *
%matplotlib inline
from IPython.display import Markdown, display

```

**Bertand**

First we wanna solve the case of Bertand competetion, i.e a sitatuion, where the firms compete in prices. Here we consider two companies who are selling differentiated products. The quantity each firm produces is dependent on the on it's own price $p_1$ and the competitor's price $p_2$.
$$x_1 = a - p_1 +bp_2$$
$$x_2 = a - p_2 + bp_1$$
The individual companies profits can be defined as: 
$$\pi_1 = (p_1 - c_1)(a - p_1 +bp_2)$$
$$\pi_2 = (p_2 - c_2)(a - p_2 + bp_1).$$
For the following parts we assume that both companies have identical marginal costs and thus $$c_1 = c_2 = c.$$


```python
#Variables
a = sm.symbols('a')
p1 = sm.symbols('p1') 
p2 = sm.symbols('p2') 
c = sm.symbols('c') 
b = sm.symbols('b') 
alpha=sm.symbols('alpha')

#Frist we define the profit funtioncs for the two firms:  
pi1 = (p1-c)*(a-p1+b*p2)
pi2 = (p2-c)*(a-p2+b*p1)
```

In order to get the equilibrium price we first have to determine the first order conditions (FOC) by maximizing the profit functions of the individual firms with repsect to the prices they offer. 
$$\frac{\partial\pi_1}{\partial p_1}= a - 2p_1 + bp_2 + c=0$$
$$\frac{\partial\pi_2}{\partial p_2}= a - 2p_2 + bp_1 + c=0$$

By solving the first order conditions for the respective prices we arrive at the reaction functions. 
$$p_1^{*} = \frac{a+ b*p_2 + c}{2}$$
$$p_2^{*} = \frac{a+ b*p_1 + c}{2}$$

By setting them equal, since it is the best repsonse for both firms to set their price equal to the competitior's price, we arrive at our equilibrium price: 
$$p^{*} = \frac{-(a+c)}{b-2}$$


```python
#The first order conditions are then: 
FOC1 = sm.diff(pi1, p1)
FOC2 = sm.diff(pi2, p2)


#Reactionsfuncions 
RF1 = sm.solve(FOC1, p1)[0]
RF2 = sm.solve(FOC2, p2)[0]

#Equilibrium price
RF12 = sm.solve(FOC2,p1)[0] 
EP1 = sm.solve(RF1-RF12,p2)[0]
print(EP1) #Equilibrium price both companies will charge 
RF21 = sm.solve(FOC1,p2)[0] 
EP2 = sm.solve(RF2-RF21,p1)[0]
```

    -(a + c)/(b - 2)


Using the lambdify feature from the sympy package we can now make a function to solve the bertrand competition for chosen values. In this case we set $a$ equal to 101, $b$ equal to 0.1 and $c$ equal to 1 (the numbers were decided on arbitrarily). This yields an equilibrium price of 53. 



```python
#make it a function to solve it 
bertrand = sm.lambdify((a,b,c), EP1)
price1=bertrand(101,0.1,1)
quantity1 = 101 - 0.9*bertrand(101,0.1,1)
profit1 = (bertrand(101,0.1,1)-1)*(101 - 0.9*bertrand(101,0.1,1))

#Defining a way to print the solution in a fancy way
def fancy(string):
    display(Markdown(string))


fancy(f'The equilibrium price is $p$ = {price1:.0f}, the equilibrium quantity is $q$ = {quantity1:.0f} and the profits for both firms are $π$ = {profit1:.0f}') 
```


The equilibrium price is $p$ = 54, the equilibrium quantity is $q$ = 53 and the profits for both firms are $π$ = 2776


**Stackelberg**

We now turn to the Stackelberg, i.e a sequantiacal competion in quantities, thus one company is the leader and the other is a follower. 


```python
#values 
p = sm.symbols('p')
x1 = sm.symbols('x1')
x2 = sm.symbols('x2')

#Giving values to our parameters
a = 101
c = 1
b = 0.1
```

Each firm chooses a quantitiy $x_{1}, x_{2}\geq0$ where the cost of production is $c_{1}=c_{2}$. Each firm set the price $p = a-b\left(x_{1}+x_{2}\right)$. The profit function for firm $i$ is then written as $\pi_{i}\left(x_{i}\right)=p\left(x_{i}\right)-c_{i}x_{i}$ Assume firm 1 is the leader and 2 is the follower. The game is solved by backward induction.
  
We will first define the inverse demand function $p$:


```python
#Defining price
p = a - b*(x1+x2)  
```

Then, we set up the profit function for firm 1 and firm 2:


```python
#Profit function of firm 1
profit_1 = p*x1 - c*x1
profit_2 = p*x2 - c*x2
```

We will say that firm 1 decides first on its output. Firm 2 observes the decision of firm 1 and decides how to reacht. We can solve the problem backwards and find optimal output of firm 2 for each possible output level of firm 1.  
It means that we maximize the profit of firm 2 by differentiating it by $x_2$.


```python
#Differentiating the profit of firm 2 by q2

diff_profit_2 = sm.diff(profit_2, x2)
diff_profit_2
```




$\displaystyle - 0.1 x_{1} - 0.2 x_{2} + 100$




```python
#Finding optimal q2 by isolating q2 from the differentiation
opt_2 = sm.solve(diff_profit_2, x2)
opt_2[0]
```




$\displaystyle 500.0 - 0.5 x_{1}$



Then, we can substitute $q_2$ into the profit function of firm 1. Thanks to that, we have a single variable function.


```python
#Substituting q2 into the profit function of firm 1
profit_1_subs = profit_1.subs(x2,opt_2[0])
profit_1_subs
```




$\displaystyle x_{1} \left(51.0 - 0.05 x_{1}\right) - x_{1}$



Once again, we need to differentiate the profit of firm 1 by $q_1$ to find the optimal value of $q_1$.


```python
#Differentiating the profit of firm 1 by q1

diff_profit_1 = sm.diff(profit_1_subs, x1)
diff_profit_1
```




$\displaystyle 50.0 - 0.1 x_{1}$




```python
#Finding the optimal value of q1

opt_1 = sm.solve(diff_profit_1, x1)
opt_1[0]
```




$\displaystyle 500.0$



Now, we can find the optimal value of $q_2$.


```python
opt_2_subs = opt_2[0].subs(x1,opt_1[0])
opt_2_subs
```




$\displaystyle 250.0$



We can also find the equilibrium market price, profit of firm 1 and profit of firm 2:


```python
#Finding equilibrium market price
p = a - b*(opt_1[0] + opt_2_subs)

#Finding profit of firm 1
profit_1 = p*opt_1[0] - c*opt_1[0]

#Finding profit of firm 2
profit_2 = p*opt_2_subs - c*opt_2_subs


fancy(f'The equilibrium price is $p$ = {int(p)}, and the profits for both firms are $π_1$ = {int(profit_1)} and $π_2$ = {int(profit_2)}')
```


The equilibrium price is $p$ = 26, and the profits for both firms are $π_1$ = 12500 and $π_2$ = 6250


We now want to make a figure that show how the beste responses relate.


```python

```


```python
#Defining the x-axis
qu_1 = np.linspace(0, 1000, 1000)

#Defining the lambda functions for our best responses functions
BR2 = lambda qu_1 : 500 - qu_1/2
BR1 = lambda qu_1 : 500 - qu_1/2

#Plotting the figure
plt.figure(figsize=(15,10))
plt.plot(qu_1, BR2(qu_1), label = "Best Response of Firm 1")
plt.plot(BR1(qu_1),qu_1 , label = "Best Response of Firm 2")
plt.grid(True)         
plt.xlabel('$q1$')
plt.ylabel('$q2$')
plt.title('Best response functions')
plt.legend()

#To include a quantity of 0 for 1 firm, we set the axis like this
plt.xlim(0,1000) # sets the x-axis
plt.ylim(0,600) # Sets the y-axis


plt.show()
```

This tells that when firm chooses a quantity of 500 in the first stage, firm 2 will choose a quantity half the size of firm 1 in second stage. In a case where firm 1 chooses to produce 1000 units, firm 2 will then produce 0. 

**Cournot**

For the Cournot case we are gonna set up some functions to solve the model in order to make an extention with n firms later on.

Consider the linear invese demand curve $p_{i}\left(x_{i},x_{j}\right)=a-x_{i}-bx_{j}$ where $b\in\left(0,1\right)$ is the elasticity of substitution between the goods, i.e. indicates that the goods are heterogenous for $0<b<1$. 

We write the inverse demand:


```python
def price(x1,x2,a,b):
    """A function that computes prices 
    args:
    x1 : quantity for firm 1 
    x2 : quantity for firm 2 
    a :total demand 
    b : elasticity of substitution
    """
    return a-b*(x1+x2)

```

Firms have constant marginal costs, where total costs are given by $c\left(x_{i,j}\right)=cx_{i,j}$


```python
def cost(x,c):
    """A function that computes total cost and sets it equal to 0 if nothing is produced
    args:
    c = cost for firm i
    x = quantity for firm i
    """
    if x == 0:
     cost = 0
    else:
     cost = c*x
    return cost
```

Firms maximize profits given the choice of the other firm such that $$\pi\left(x_{i},x_{j}\right)=p_{i}\left(x_{i},x_{j}\right)x_{i}-c\left(x_{i}\right)$$


```python
def profit(x1,x2,c1,a,b):
    """A function that computes profit for firm 1
    args:
    x1 : quantity for firm 1 
    x2 : quantity for firm 2 
    c1 : cost for firm 1 
    a : total demand 
    b : elasticity of substitution
    """
    return price(x1,x2,a,b)*x1-cost(x1,c1)
```

In a Cournot setting, firm $i$ chooses $x_{i}$ taking $x_{j}$ as given. Thus, the Nash equilibrium for $x_{i}^{*},x_{j}^{*}$ holds for $$x_{i}^{*}=\text{arg}\underset{x_{i}}{\text{ max }}\pi\left(x_{i},x_{j}^{*}\right)\;\;\;\;\;\;\text{for each }i\neq j$$ Thus, we are able to find firm $i$ optimal production level given the choice of $x_{j}$ by creating the function 'reaction'. Hence, taking the first order condition to find the Nash equilibrium for $x_{i}^{*},x_{j}^{*}$ yields $$\frac{\partial\pi\left(x_{i},x_{j}\right)}{\partial x_{i,j}}=0$$

In this case, we also want to calculate the profit maximizing outcome, so we use the function 'brute'. This allows one to look for the global maximum of $\pi$. Likewise, we keep in mind that 'scipy optimize' is bound to minimize problems, and since we want to maximize the profit function, we then minimize $-\pi\left(x_{i},x_{j}\right)$. 


```python
def reaction(x2,c1,b,a):
    """A function that computes the optimal quantity
    args
    x2 : quantity for firm 2
    c1 : cost for firm 1
    b : elasticity of substitution
    a : total demand 
    """
    
    #bounds = ((0,a/b-x2),)
    x1 = optimize.minimize(lambda x: -profit(x,x2,c1,a,b),x0 = [x2],method='L-BFGS-B')
    #print(x1.message)
                                                                     
    return x1.x[0]
```


```python
def reaction(x2,c1,b,a):
    """A function that computes the optimal quantity
    args
    x2 : quantity for firm 2
    c1 : cost for firm 1
    b : elasticity of substitution
    a : total demand 
    """
    x1 = optimize.brute(lambda x: -profit(x,x2,c1,a,b), ((0,1,),)) 
                                                                 
    return x1[0]
```


```python
reaction(333.3715739255503,1,0.1,101)
```




    333.31425781250056



Now, we must find an equilibrium for price and quantity given their best response. For this reason, consider the vector function $f\left(x_{i},x_{j}\right)$ which takes firm $i$'s best response $r_{i}\left(x_{j}\right)$ to $j$'s production level $x_{j}$ into consideration. We can write this as $$f\left(x_{i},x_{j}\right)=\left[\begin{array}{c}
r_{i}\left(x_{j}^{*}\right)\\
r_{j}\left(x_{i}^{*}\right)
\end{array}\right]$$ 

We are then looking for the quantities produced is equal to best response given what the other firm produce, i.e the point where: 

$$\left(\begin{array}{c}
x_{1}^{*}\\
x_{2}^{*}
\end{array}\right)=\left(\begin{array}{c}
r_{1}(x_{2}^{*})\\
r_{2}(x_{1}^{*})
\end{array}\right)$$

Thus, we create a function 'vec_reaction' given by the difference between an array of the quantaites subtracted by the best response functions and pass $b,a, c_{i},c_{j}$ as a vector of parameters. 


```python
def vec_reaction(x,param): 
    """A function that computes the difference between the quantities and the beste response functions
    args
    x : quantaties
    param : list of parameters"""
    return array(x)-array([reaction(x[1],param[2],param[0],param[1]),reaction(x[0],param[3],param[0], param[1])])
```

Now, moving on to find the actual Cournot equilibrium, we are able to set the values of each paramter. Using the 'fsolve' from 'scipy.optimize' to find $x$. Also, giving 'fsolve' the inital guess of the quantaites, while passing the paramater vector 'param' as arguments 'args' to 'fsolve', this yields:


```python
param = [0.1,101,1,1]
x0 = [100,100]

EQ = optimize.fsolve(vec_reaction, x0, args = (param))
EQ

fancy(f'The Nash equilibrium quantity for both firms are $x^* =[x_1^*,x_2^*]$ = {EQ} where the equilibrium price for both firms are $p^*$ = {price(EQ[0],EQ[1],101,0.1)} with a profit of $\pi^*$= {profit(EQ[0],EQ[1],1,101,0.1)}')

```


The Nash equilibrium quantity for both firms are $x^* =[x_1^*,x_2^*]$ = [333.33330078 333.33330078] where the equilibrium price for both firms are $p^*$ = 34.333339843749854 with a profit of $\pi^*$= 11111.11219618032


We now want to show, how the quantity produced changes, when the cost of firm 1 changes.


```python
range_c = np.linspace(1,0,20)
range_x = [optimize.fsolve(vec_reaction, x0, args = ([0.1,101, c, 1])) for c in range_c]

plt.style.use('seaborn')
plt.clf()
plt.plot(range_c,array(range_x)[:,0], label = '$x_1$')
plt.plot(range_c,array(range_x)[:,1], label = '$x_2$')
plt.xlabel('$c_1$')
plt.ylabel('$x$')
plt.legend()
plt.show()
```

In the figure, we see, that as the cost of firm 1 goes up, the quantity produced by firm 1 falls, while the quantity produced by firm 2 goes up.

**Compare all three cases**

We have now tried to solve a situation with two firms competing with three different models, giving them the same parameter values. This leads to the following results:


```python
fancy(f'For Bertrand: The equilibrium price is $p$ = {int(price1)}, the equilibrium quantity is $q$ = {int(quantity1)} and the profits for both firms are $π$ = {int(profit1)}') 
fancy(f'For Stackelberg: the equilibrium price is $p$ = {int(p)}, and the profits for both firms are $π_1$ = {int(profit_1)} and $π_2$ = {int(profit_2)}')
fancy(f'For Cournot: the equilibrium quantity for both firms are $x^* =[x_1^*,x_2^*]$ = {EQ} where the equilibrium price for both firms are $p^*$ = {price(EQ[0],EQ[1],101,0.1)} with a profit of $\pi^*$= {profit(EQ[0],EQ[1],1,101,0.1)}')

```


For Bertrand: The equilibrium price is $p$ = 53, the equilibrium quantity is $q$ = 52 and the profits for both firms are $π$ = 2775



For Stackelberg: the equilibrium price is $p$ = 26, and the profits for both firms are $π_1$ = 12500 and $π_2$ = 6250



For Cournot: the equilibrium quantity for both firms are $x^* =[x_1^*,x_2^*]$ = [333.33330078 333.33330078] where the equilibrium price for both firms are $p^*$ = 34.333339843749854 with a profit of $\pi^*$= 11111.11219618032


**Cournot with n> 2 firms**

In order to do an extention to one of the models, we now want to investigate the case of N > 2 firms. First we set up the situation theoretically. In the following theoretical set-up we will follow the set-up in this note, adapting it to our situation: http://people.exeter.ac.uk/dgbalken/micro07/CournotHandout.pdf (we add the posibility of different cost for each firms)

Following the notation from ealier we now add: 

- Firm $i$'s output: $x_i$

- Total output in in industry: $X=x_1+x_2+...+x_n$

- The opponent's output: $x_{-i} = X-x_i$



Each firm will then max profit given the expectation of $x_{-i}$: 

$$\Pi_i(x_{-i},x_i)=(p(x_{-i}+x_i)-c_i)x_i$$

as with the n = 2 case we use the invese demand funcktion in a linear case:

$$p = a-b(x_{-i}+x_i)=a-bX$$ 

The FOC of is then:
$$\frac{\partial \Pi_i}{\partial x_i} = \frac{\partial p}{\partial x_i} \cdot x_i + p-c_i = 0$$


where $\frac{\partial p}{\partial x_i} = -b$ giving us

$$\frac{\partial\Pi_i}{\partial x_i} : -bx_i + p-c_i = 0$$

inserting the inverse demand function:

$$\frac{\partial\Pi_i}{\partial x_i} : -bx_i + (a-b(x_{-i}+x_{i}))-c_i = 0$$

solving for $x_i$ to get the reaction function for firm $i$ given its expectation for the amount the rest of the firms produce:

$$x_i(x_{-i}) = \frac{a-c}{2b}-\frac{1}{2}x_{-i}$$

which is the case for all firms $i=1,..,n$. We have $X=x_i+x_{-i}$ so the FOC found earlier can be rewritten to:
$$ -bx_i + (a-bX-c_i) = 0$$

which for the n firms:
$$-bX + n(a-bX)-n\bar c = 0$$

where we use the average cost level: $\bar c = \frac{c_1+c_2+...+c_n}{n}$

we can now rewrite this and get the quantity produced and the price in the market: 

$$bX + n(bX) = n(A-\bar c)$$

$$(n+1)bX  = n(A-\bar c)$$

Quantity is then:

$$X = \frac{n}{n+1} \frac{a-\bar c}{b}$$

Price is:
$$p = a-bX = a-b[\frac{n}{n+1} \frac{a-\bar c}{b}] = \frac{1}{n+1}a + \frac{n}{n+1}\bar c$$


Using the found functions we can now try to set a program that is able to solve cases with N firms. 



```python
def Round(a): 
    """A function that prints floats rounded as strings
    args:
    a : a float"""
    return str("{:.0f}".format(a))
```

First we define a function, that calculates the prices $p = \frac{1}{n+1}a + \frac{n}{n+1}\bar c $ rewriting using the definition of the average cost level gives: $p= \frac{1}{n+1}a + \frac{c_1+c_2+...+c_n}{n+1}$ which we use in our function.


```python
def price_eq(a, n, c): #Returns the equilibirum price given n firms 
    """A function that computes price
    args:
    a : total demand 
    n : number of firms 
    c : cost of firm i 
    """
        
    price = max((a+sum(c))/(1+n), 0)
    return price
```

We now want to define a function that returns the quantity produced for firm $i$.
this is given by:
$$x_i = \frac{a-bX-c_i}{b}=\frac{a-c_i}{b}-\frac{n}{n+1}\frac{a-\bar c}{b} = \frac{1}{n+1}\frac{a}{b}+ \frac{n(\bar c -c_i)-c_i}{(n+1)b} = \frac{1}{n+1}\frac{1}{b}[a- n(\bar c -c_i)-c_i]$$

use that $\bar c = \frac{c_1+c_2+...+c_n}{n}$
$$x_i = \frac{1}{n+1}\frac{1}{b}[a+ nc_i+(c_1+c_2+...+c_n-c_i)]$$
which we now use to define a function, that returns the quantity: 


```python
def xi_eq(a, b, n, c, i): #Returns equilibrium quantity for firm i
    """A function that computes quantity for firm i
    args:
    a (float):total demand 
    b (float):elasticity of substitution between the good
    n (float):number of firms 
    c (float): cost of firm i 
    i (integer): the number of the firm
    Returns a float"""
    
    #x_i = (1/b)*((a+n*c[i]+(sum(c)-c[i]))/(1+n))
    #if x_i < 1: 
     #   x_i = 0
    x_i = (1/b)*((a+(n-1)*c[i]-sum(c))/(1+n))   
    return x_i
```

We can then define a function, that returns the quantity for the entire number of firms 


```python
def X_eq(a, b, n, c): 
    """A function that computes quantity for the total amount of firms 
    args:
    a : total demand 
    b : elasticity of substitution between the good
    n : number of firms 
    c : cost of firm i 
    """
    
    X = []
    for i in range(0, n, 1):
        X.append(xi_eq(a, b, n, c, i))
        
    return sum(X)
```


```python

```


```python
n= 10
c_vec = [3 for i in range(n)]
X_eq(a,b,n,c_vec)
```




    890.909090909091




```python
len(c_vec)
```




    10




```python
n= 1000
c_vec = [random.uniform(1,3) for i in range(n)]
X_eq(a,b,n,c_vec)
```




    989.1567070866106




```python
np.sum([a/b/(n+1)+(n*(np.mean(c_vec)-c)-c)/((n+1)*b) for c in c_vec])
```




    989.1567070866333




```python
n/(n+1) * (a-np.mean(c_vec))/b
```




    989.1567070866319



Based on the assumtion that $$c_i = \bar c = c $$ the profit for firm $i$ can be found by:
$$\Pi_{i}=(p-c)x_{i}=\left(\frac{a+nc}{n+1}-c\right)\frac{1}{n+1}\frac{1}{b}(a-c)$$

Without identical marginal costs $$c_i \neq c $$ $$\bar c \neq c$$ the profit for firm $i$ would be: 
$$\Pi_{i}=(p-c_i)x_{i}=\left(\frac{a+n\bar c}{n+1}-c_i\right)\frac{1}{n+1}\frac{1}{b}(a+n(\bar c - c_i)-c_1) = \frac{1}{(n+1)^2 * b}(a+n(\bar c - c_1)-c_i)^2$$

By maximizing the previously derived profit function we find the equilibrium profit.


```python
def profit_eq(a, b, n, c, i): 
    """A function that computes profit for firm i and sets profit to zero, if there is no production
    args:
    a : total demand 
    b : elasticity of substitution between the good
    n : number of firms 
    c : cost of firm i 
    i : the number of the firm
    """
    
    if xi_eq(a, b, n, c, i) == 0: 
        profit = 0
    else:
        profit = max((1/b)*(((a-n*c[i]+(sum(c)-c[i]))/(1+n))**2), 0)
        
    return profit

```

We are now ready to define a function, that can solve a general case of cournot competetion with n firms. We choose the number of firms to be drawn randomely using random.randint, while the cost for each firm will be drawn from a uniform destribution using random.uniform. 


```python
def EQ_N(a, b, N, c_low, c_high): 
    """A function that draws number of firms and compute quantaties, price and average profit
    args:
    a (float):total demand 
    b (float):elasticity of substitution between the good
    N (float): the maximal number of firms that can be drawn 
    c_low (float): the lower bound of the interval for the cost 
    c_high (float): the upper bound of the interval for the cost 
    Returns a float"""
 
  
    #draw the number of firms:
    n = random.randint(2, N) 
    
    #drawing of costs for each firm:
    c = []
    for i in range(0, n+1, 1):
        c_i = random.uniform(c_low, c_high)
        c.append(c_i)
    
    #the price: 
    price = price_eq(a, n, c) 
    
    #quantities and profits
    xi_list = []
    profiti_list = []
    for i in range(0, n+1, 1):
        xi_list.append(xi_eq(a, b, n, c, i))
        profiti_list.append(profit_eq(a, b, n, c, i))
   
    
    #the quantity for the entire number of firms
    X = X_eq(a, b, n, c)

    
        
    output = print('\n The number of firms drawn are ' + Round(n)  + 
                  '\n The total amount produced is ' + Round(X) +',' + ' while the price is ' + Round(price) + 
                   '\n The average profit for each firm is ' + Round(np.array(profiti_list).mean()) + '.'
                   '\n')
    
    return output
    

    
```

We are now ready to try a case with parameters as in the n = 2, but with a randomw number of firms drawn:


```python
a = 101
b = 0.1
N = 25
c_low = 1
c_high = 5
EQ_N(a, b, N, c_low, c_high)

```

    
     The number of firms drawn are 9
     The total amount produced is 848, while the price is 13
     The average profit for each firm is 1032.
    



```python
def EQ_N(a, b, n, c_low, c_high): 
    """A function that draws number of firms and compute quantaties, price and average profit
    args:
    a (float):total demand 
    b (float):elasticity of substitution between the good
    N (float): the maximal number of firms that can be drawn 
    c_low (float): the lower bound of the interval for the cost 
    c_high (float): the upper bound of the interval for the cost 
    Returns a float"""
 
  
    #draw the number of firms:
    #n = random.randint(2, N) 
    
    #drawing of costs for each firm:
    c = []
    for i in range(0, n, 1):
        c_i = random.uniform(c_low, c_high)
        c.append(c_i)
    
    #the price: 
    price = price_eq(a, n, c) 
    
    #quantities and profits
    xi_list = []
    profiti_list = []
    #for i in range(0, n+1, 1):
        #xi_list.append(xi_eq(a, b, n, c, i))
        #profiti_list.append(profit_eq(a, b, n, c, i))
   
    
    #the quantity for the entire number of firms
    X = X_eq(a, b, n, c)

    return X,price, c
    output = print('\n The number of firms drawn are ' + Round(n)  + 
                  '\n The total amount produced is ' + Round(X) +',' + ' while the price is ' + Round(price) + 
                   '\n The average profit for each firm is ' + Round(np.array(profiti_list).mean()) + '.'
                   '\n')
    
    return output
    

    
```


```python
X,price, c_vec =  EQ_N(a, b, 1000, c_low, c_high)
print(X,price)
```

    979.4548976585473 3.054510234142361



```python
n_vec = list(range(100,1000))


output_vec1 = array([EQ_N(a,b,n,c_low,c_high) for n in n_vec])


fig,ax = plt.subplots()
ax.plot(n_vec,output_vec1[:,0], label = '$X$')
plt.legend()
plt.show()
```


```python

fig,ax = plt.subplots()
ax.plot(n_vec,output_vec1[:,1], label = '$X$')
plt.legend()
plt.show()
```


```python
def solution(a,b,c_vec):
    n= len(c_vec)
    q = n/(n+1) * (a-np.mean(c_vec))/b
    p = a-b*q
    return q,p
```


```python
solution(a,b,c_vec)
```




    (979.4548976585766, 3.0545102341423416)




```python
solution(a,b,[1,2,3])
```




    (742.5, 26.75)




```python
def setup_solve(a,b,n,c_low,c_high):
    c_vec = [random.uniform(c_low, c_high) for i in range(n)]
    return solution(a,b,c_vec)
```


```python
n_vec = list(range(100,1000))
c = 3

#output_vec = array([setup_solve(a,b,n,1,5) for n in n_vec])
output_vec = array([solution(a,b,c_vec) for c_vec in output_vec1[:,2]])


fig,ax = plt.subplots()
ax.plot(n_vec,output_vec[:,0], label = '$X$')
plt.legend()
plt.show()
```


```python
plt.style.use('seaborn')
plt.clf()
plt.plot(n_vec,output_vec[:,1], label = '$price$')
plt.legend()
plt.show()
```


```python
output_vec[-1,1]
```




    3.129127017124361



**Solving using python**



```python
N = 100
c_low, c_high = 3,5
c_vec = np.random.uniform(c_low,c_high,size=N)
```


```python
def focs(a,b,q_vec,c_vec):
    # Using analytical FOCs from pdf-link
    return -b*q_vec+a-b*np.sum(q_vec)-c_vec
```


```python
def solve(a,b,c_vec):
    obj = lambda q_vec : focs(a,b,q_vec,c_vec)

    res = optimize.root(obj, x0 = [90 for n in range(N)])
     
    q_vec_star = res.x
    Q = np.sum(q_vec_star)
    
    # Use inverse demand to find price 
    p = a- b*Q
    return Q, p
    
    
```


```python
solve(a,b,c_vec)
```




    (960.2800162703653, 4.971998372963469)




```python
# Compare
solution(a,b,c_vec)
```




    (960.2800162703655, 4.971998372963441)



**Conclusion**

In this paper we examine two firms who have the possibility of competing in Bertrand, Stackelberg and Cournot. We use the same parameter values for all competition types, and find that the firms benefit the most from competing in a Cournot setting as this yields the highest total profits of 22,222. This is intuitive, as in contrast to Bertrand, the firms will undercut each other until their price is equal to their marginal costs to obtain efficiency. Therefore, this setting yields a price as in perfect competition, whereas Cournot yields an oligopoly equilibrium. 

As Cournot yields the best outcome in terms of profits for two firms, we extend this setting for $n$ firms. We find that for $n\rightarrow\infty$, the average profit for a firm decreases because competing in quantities becomes more intense as the competition setting moves from oligopoly to infinitely many firms competing in quantities.

c2bb25b87b895dfebeda7f5c4941520c4e268180
ipynb
Jupyter Notebook
modelproject/modelproject_done4.0.ipynb
AskerNC/projects-2021-aristochats
cade4c02de648f4cd1220216598dc24b67bb8559
MIT
modelproject/modelproject_done4.0.ipynb
AskerNC/projects-2021-aristochats
cade4c02de648f4cd1220216598dc24b67bb8559
MIT
modelproject/modelproject_done4.0.ipynb
AskerNC/projects-2021-aristochats
cade4c02de648f4cd1220216598dc24b67bb8559
MIT
Qwen/Qwen-72B
1. YES
2. YES