How to reduce the model generation time?

Hello,

The generation time for the below model consisting of 112629 binary variables and 398918 constraints is more than 2500 seconds. It is longer than I would expect. I am wondering why it takes this much time to generate the model. The Python code can be found below. Is there a problem related to my code? How can I improve the efficiency? 

Thank you for your ideas in advance.

import gurobipy as gp
from gurobipy import GRB
import random
import numpy as np
from datetime import datetime
rand=random.Random()
def createList(r1, r2): 
    return [item for item in range(r1, r2+1)]
import pandas as pd
rand.seed(2353573)
#Define and initialize parameters related to order_type 
sku_num=328
skus=createList(0, sku_num-1)
additional_sku_num=300
additional_skus=createList(0, additional_sku_num-1)
total_sku_num=sku_num+additional_sku_num
total_skus=createList(0, total_sku_num-1)
#Read the file order_type.csv including a binary matrix with 246912 rows and 328 columns
chunk= pd.read_csv('order_type.csv',header=None,chunksize=10000)
order_type2 = pd.concat(chunk)
order_type2= order_type2[0].str.split(";",expand=True,)
order_type2[skus] = order_type2[skus].astype(int)
#Read the file prob.csv including 1 row and 328 columns
prob = pd.read_csv('prob.csv',header=None)
prob= prob[0].str.split(";",expand=True,)
prob[skus] = prob[skus].astype(float)
#Generate new random orders
additional_orders = pd.DataFrame(0, index=np.arange(len(order_type2)), columns=additional_skus)
additional_prob = pd.DataFrame(0, index=np.arange(len(prob)), columns=additional_skus)
additional_prob[additional_skus] = additional_prob[additional_skus].astype(float)
for s in range(additional_sku_num):
    additional_prob.iloc[0,s]=rand.uniform(0,1)
for o in range(len(additional_orders)):
    if rand.uniform(0,1)<=0.75:
        degree_to_be=1
    elif rand.uniform(0,1)>0.75 and rand.uniform(0,1)<=0.9:
        degree_to_be=2
    elif rand.uniform(0,1)>0.9 and rand.uniform(0,1)<=0.97:
        degree_to_be=3
    else:
        degree_to_be=4
    cttot=0
    while cttot<degree_to_be:
        considered_sku=rand.randint(0,additional_sku_num-1)
        if additional_orders.iloc[o,considered_sku]==0 and rand.uniform(0,1)<=additional_prob.iloc[0,considered_sku]:
            additional_orders.iloc[o,considered_sku]=1
            cttot=cttot+1            
#Concatenate given orders and the additional ones. Count and delete duplicates
order_type3 = pd.concat([order_type2,additional_orders], axis=1, ignore_index=True)
order_type=order_type3.groupby(order_type3.columns.tolist(),as_index=False).size()
amount=order_type["size"]
del order_type["size"] #Size of the order_type=(112001,628)
degree = np.sum(order_type,axis=1).tolist()
order_num=len(order_type) #order_num=246912
orders=createList(0, order_num-1)
capacity=62
now = datetime.now()
jehl = gp.Model('jehl')
x=jehl.addVars(total_sku_num, vtype=GRB.BINARY, name="x")
y=jehl.addVars(order_num, vtype=GRB.BINARY, name="y")
con1 = jehl.addConstr((gp.quicksum(x[s] for s in total_skus) <= capacity))
con2 = jehl.addConstrs((y[o] <= x[s] for o in orders for s in total_skus if order_type.iloc[o,s]==1))
jehl.setObjective(( gp.quicksum(amount.iloc[o]*y[o] for o in orders)),
                GRB.MAXIMIZE)
jehl.setParam('LogToConsole', 0)
jehl.setParam('TimeLimit',3600)
jehl.setParam('MIPgap',0)
generation_time=(datetime.now()-now).total_seconds()
jehl.optimize()