Finance

import numpy as np
from cvxopt import matrix, solvers
import matplotlib.pyplot as plt
import math 


n=11
m=21
a=-0.5
b=0.7
p=b/(b-a)

from cvxopt import printing
printing.options["width"]=-1
printing.options["height"]=-1
    
def cout(x,y):
    if y<x:
        return 0
    return (y-x)


def creerA(n,m):
    A=[]
    for i in range(n):
        ligne=[0. for k in range(n*m)]
        for j in range(m):
            ligne[i*m+j]=1.
        A.append(ligne)
    for j in range(m-1):
        ligne=[0. for k in range(n*m)]
        for i in range(n):
            ligne[i*m+j]=1.
        A.append(ligne)
    A=np.array(A)
    return A
    
def creerA2(n,m):
    A=[]
    for i in range(n-1):
        ligne=[0. for k in range(n*m)]
        for j in range(m):
            ligne[i*m+j]=(1+a)**j*(1+b)**(m-1-j)
        A.append(ligne)
    A=np.array(A)
    return A
    
def creerABS(A1,A2):
    A=np.concatenate((A1,A2),axis=0)
    return A
    


def solve(c,d,n,m):    
    A1=matrix(creerA(n,m))
    A2=matrix(creerA2(n,m))
    A=matrix(creerABS(A1,A2))
    # print(A.size)
    c=matrix(c)
    # print(c)
    d=matrix(d)
    # print(d)
    # print(d.size)
    h=np.zeros((n*m))
    h=matrix(h)
    G=-np.eye((n*m))
    G=matrix(G)
    # print(A)
    # print(G)
    dims ={'l':G.size[0], 'q':[], 's':[]}
    return solvers.conelp(c,G,h,dims,A,d)

def reassemblage(PI,n,m):
    M=np.zeros((n,m))
    for i in range(m):
        for j in range(n):
            M[j,i]=PI[i+j*m]
    return M


d=np.zeros((2*n+m-2))

for i in range(n):
    d[i]=(math.factorial(n-1) / (math.factorial(i)*(math.factorial(n-1-i))))*((1-p)**(n-1-i))*p**(i)
    
    
for j in range (m-1):
    d[n+j]=math.factorial(m-1) / (math.factorial(j)*(math.factorial(m-1-j)))*((1-p)**(m-1-j))*p**(j)
    
for i in range(n-1): 
    d[n+m-1+i]=(math.factorial(n-1) / (math.factorial(i)*(math.factorial(n-1-i))))*(1-p)**(n-1-i)*p**(i)*((1+a)**i)*((1+b)**(n-1-i))
    
c=np.zeros(n*m)
        
for i in range(n):
    for j in range(m):
        c[j+m*i]=cout(((1+a)**i)*(1+b)**(n-1-i),((1+a)**j)*(1+b)**(m-1-j))
        
PI=solve(c,d,n,m)['x']
min=0
for i in range (n*m):
    min=min+PI[i]*c[i]
    

M=reassemblage(PI,n,m)
plt.imshow(M)

for i in range (n*m):
    c[i]=-c[i]
PI2=solve(c,d,n,m)['x']
max=0
for i in range (n*m):
    max=max+PI2[i]*c[i]
    
print('min=',min)
print('max=', -max)

# Vérification de l'accord avec le modèle CRR

x=np.zeros(n)
px=np.zeros(n)
y=np.zeros(n*m)
py=np.zeros(n*m)

for i in range(n):
    x[i]=((1+a)**i)*((1+b)**(n-i-1))
    px[i]=(math.factorial(n-1) / (math.factorial(i)*(math.factorial(n-1-i))))*((1-p)**(n-1-i))*p**(i)

for i in range(n):
    for j in range(m):
        y[i*m+j]=cout(((1+a)**(j)*((1+b)**(m-j-1))),x[i])
        if i>j:
            py[i*m+j]=0
        else :
            if (j-i)>(m-n):
                py[i*m+j]=0
            else:
                py[i*m+j]=(math.factorial(m-n) / (math.factorial(j-i)*(math.factorial(m-n-j+i))))*(math.factorial(n-1) / (math.factorial(i)*(math.factorial(n-1-i))))*((1-p)**(m-1-j))*(p**(j))
 
 
s=0
for i in range(n):
    for j in range(m):
        s=s+y[i*m+j]*py[i*m+j]
        
        
print("le coût de l'option dans le modèle CRR est",s)