BYU-PRISM / GEKKO

GEKKO Python for Machine Learning and Dynamic Optimization
https://machinelearning.byu.edu
Other
573 stars 102 forks source link

Error: At line 463 of file custom.f90: Memory allocation failed #133

Closed JustinNeumann closed 2 years ago

JustinNeumann commented 2 years ago

Dear maintainer,

I have the pleasure to use your library recently, and am facing an issue when scaling up my optimisation model. Here is a minimal example which results in the stated error. You can adjust the number of entities to be processed at the top of the code snippet.

I appreciate the help! Thank you a lot :-)

import numpy as np
import pandas as pd
import datetime
from tqdm import tqdm
from scipy.stats import gaussian_kde
import math
import os
import matplotlib.pyplot as plt
from scipy.stats import norm

from gekko import GEKKO

NUM_ENTITIES = 10 # 10 crashes already... just a few work out fine

# create sample data
di_kde = dict()
for entity in range(NUM_ENTITIES):
    arr = np.rint(np.random.normal(50, 25, size=200))
    kernel = gaussian_kde(arr, bw_method= 'scott')
    di_kde[entity] = (kernel, int(np.amax(arr)), int(np.amin(arr)))

def cdf_gauss_nogecko(kdes, maxi):
    bw = kdes.neff**(-1./(1+4))
    res = np.divide(np.cumsum([sum([norm.pdf(i, val, bw) for val in kdes.dataset.flatten()]) \
                                   for i in np.linspace(0, maxi, 100, endpoint=True)]), \
                        max(np.cumsum([sum([norm.pdf(i, val, bw) for val in kdes.dataset.flatten()]) \
                                       for i in np.linspace(0, maxi, 100, endpoint=True)])))
    return res

def pdf_gauss_gekko(kdes, yi):
    bw = kdes.neff**(-1./(1+4))
    res_gecko = m.sum([normpdf_gekko(yi, val, bw) for val in kdes.dataset.flatten()])
    return res_gecko

def normpdf_gekko(x, mean, sd):
    var = float(sd)**2
    denom = (2*math.pi*var)**.5
    num = m.exp(-(x-float(mean))**2/(2*var))
    return num/denom

cost = 0.1
revenue = 1

print(f'{datetime.datetime.now().strftime("%H:%M:%S ")}Running optimisation under constraint...')

if not os.path.isdir(os.path.abspath(r'.\Logging')):
    os.mkdir(os.path.abspath(r'.\Logging'))

m = GEKKO(remote=False) # Initialize gekko

m._path = os.path.abspath(r'.\Logging')

m.options.SOLVER = 3  # APOPT (1) is an MINLP solver, IPOPT (3) an NLP solver
m.options.IMODE = 2 # 3 or 6
m.options.MAX_MEMORY = 8

# optional solver settings with APOPT
m.solver_options = ['minlp_maximum_iterations 500', \
                    # minlp iterations with integer solution
                    'minlp_max_iter_with_int_sol 10', \
                    # treat minlp as nlp
                    'minlp_as_nlp 0', \
                    # nlp sub-problem max iterations
                    'nlp_maximum_iterations 50', \
                    # 1 = depth first, 2 = breadth first
                    'minlp_branch_method 1', \
                    # maximum deviation from whole number
                    'minlp_integer_tol 0.05', \
                    # covergence tolerance
                    'minlp_gap_tol 0.01']

print(f'{datetime.datetime.now().strftime("%H:%M:%S ")}Calculating splines...')

max_array_length = max([v[1] for k, v in di_kde.items()])

# approximate cdf with cspline
x_nogekko = [np.linspace(0, di_kde[k][1], 100, endpoint=True) for idx1, k in enumerate(di_kde.keys())]
y_nogekko = [cdf_gauss_nogecko(di_kde[k][0], di_kde[k][1]) for idx1, k in enumerate(di_kde.keys())]
x = [m.Param(value=[j for j in range(max_array_length)]) for idx1, k in enumerate(di_kde.keys())]
y = [m.Var(value=[0 for j in range(max_array_length)]) for idx1, k in enumerate(di_kde.keys())]

# approximate cdf with cspline
[m.cspline(x[idx1], y[idx1], x_nogekko[idx1], y_nogekko[idx1], True) for idx1, k in enumerate(di_kde.keys())]

print(f'{datetime.datetime.now().strftime("%H:%M:%S ")}Adding variables...')

# Initialize variables
qi = [m.Var(lb=max(0, di_kde[k][1]-10), ub=di_kde[k][1]+10, integer=True, name=f'qi_{k}') for idx, k in enumerate(di_kde.keys())]
z = [m.Var(name=f'z_{k}') for idx, k in enumerate(di_kde.keys())]
ui = [[m.Var(lb=0, ub=1, integer=True, name=f'ui_{k}_{j}') for j in range(di_kde[k][1])] for idx, k in enumerate(di_kde.keys())]

print(f'{datetime.datetime.now().strftime("%H:%M:%S ")}Adding constraint equations...')

# add constraints
for idx1, k in enumerate(di_kde.keys()):
    m.Equation(m.sum([u for u in ui[idx1]]) == 1) # only one u is 1
    m.Equation(m.sum([u * (idx2+1) for idx2, u in enumerate(ui[idx1])]) == qi[idx1]) # only one u is 1, so qi is index
    m.Equation(m.sum([y[idx1][idx2] * u for idx2, u in enumerate(ui[idx1])] ==  z[idx1]))

# add constraint
m.Equation(m.sum([z[idx4] / len(di_kde) for idx4, _ in enumerate(di_kde.keys())]) >= 0.9014)

print(f'{datetime.datetime.now().strftime("%H:%M:%S ")}Adding maximisation function...')

for ii, kk in enumerate(di_kde.keys()):
    m.Maximize(m.sum([(-cost * (m.max3(0, qi[ii]-j) + m.max3(0, j-qi[ii]) + qi[ii]) \
                           + revenue * m.min3(qi[ii], j)) * pdf_gauss_gekko(di_kde[kk][0], j) \
                          for j in range(di_kde[kk][2], di_kde[kk][1])])) # Objective

print(f'{datetime.datetime.now().strftime("%H:%M:%S ")}Start solving...')

m.open_folder()
m.options.DIAGLEVEL=1
output = m.solve(disp=True)

print(f'{datetime.datetime.now().strftime("%H:%M:%S ")}End solving...')
print('Results')
print('qi: ' + str(qi))
print('z: ' + str(z))
print('Objective: ' + str(-m.options.objfcnval))
APMonitor commented 2 years ago

The code was missing a parenthesis. I updated the code and ran it locally on Windows. It came up with a different error with this line:

    m.Equation(m.sum([y[idx1][idx2] * u for idx2, u in enumerate(ui[idx1])] ==  z[idx1]))
Traceback (most recent call last):
  File "C:\Users\johnh\Desktop\test.py", line 99, in <module>
    m.Equation(m.sum([y[idx1][idx2] * u for idx2, u in enumerate(ui[idx1])] ==  z[idx1]))
  File "C:\Users\johnh\Python39\lib\site-packages\gekko\gekko.py", line 1540, in sum
    raise TypeError("x must be a python list of GEKKO parameters, variables, or expressions")
TypeError: x must be a python list of GEKKO parameters, variables, or expressions

Is this the error that you received when you try to run the code? There is also an active gekko community on StackOverflow. This question may reach a larger audience if posted there.

JustinNeumann commented 2 years ago

Ok, I will try to post on StackOverflow asap. Best, Justin🙋‍♂️

JustinNeumann commented 2 years ago

Here is the related link: https://stackoverflow.com/questions/70256978/error-at-line-463-of-file-custom-f90-memory-allocation-failed

APMonitor commented 2 years ago

Thanks for posting on Stack Overflow. We'll take the discussion there.