El mar, 05-06-2007 a las 09:33 -0700, J.-O. Moussafir escribió: > Hello, > > I've come across segfault using Python-glpk. Does anyone have had some experience with it ? > > Jacques-Olivier > --------------------------------------------------------------------------------------------------------------------------
> I have experience with PULP (A linear programming modeler in Python - http://www.jeannot.org/~js/code/index.en.html#PuLP). > This is a program developed in python to use GLPK in python. I send you an example of "data envelopment analysis" > developed using Pulp. ############################################################### #!/usr/bin/env python # @(#) $Jeannot: test1.py,v 1.11 2005/01/06 21:22:39 js Exp $ # Importar PuLP: A python linear programming modeler from pulp import * # Importar Scipy: Scientific python from scipy import * from string import * # Importar otras librerias import Gnuplot, Gnuplot.funcutils,os,time,datetime # Funciones a utilizar # Funcion para leer datos def LeeDat(archivo): # lee datos del archivo f = open(archivo, 'r') salida = io.read_array(f) f.close() return salida # Funcion para generar archivos def Guarda(variable,nombre): # Guarda variable en archivo archivo = file(nombre, 'w') # Creando archivo io.write_array(archivo,variable) # Escribiendo archivo archivo.close() # Cerrando archivo # Comienza el programa intro = "\n\ Data Envelopment Analysis\n\ Version Beta-0.1\n\ Carlos Herrera - Andres Vejar\n\ DEPARTAMENTO DE INGENIERIA INDUSTRIAL\n\ UNIVERSIDAD DE SANTIAGO DE CHILE\n\ ::::::::::::::::::::::::::::::::" print intro print "" print "Parametros del modelo" print "---------------------" i1 = 19 print "El numero de entidades a evaluar es",i1 print "" i2 = int(raw_input("Ingrese el numero de meses: ")) # Parametros del modelo dmu = i1 # numero de entidades meses = i2 # numero de meses # Matriz de ponderaciones de eficiencia solu0 = zeros((dmu,meses)) solu1 = mat(zeros((dmu,meses))) solu2 = mat(zeros((dmu,meses))) for mes in range(meses): entrada0 = 'entrada-' + str(mes+1) +'.datos' salida0 = 'salida-' + str(mes+1) +'.datos' entrada1 = 'entrada-' + str(mes+2) +'.datos' salida1 = 'salida-' + str(mes+2) +'.datos' datIn0 = LeeDat(entrada0) datOut0 = LeeDat(salida0) datIn1 = LeeDat(entrada1) datOut1 = LeeDat(salida1) entr = shape(datIn0)[1] sal = shape(datOut0)[1] # ANALISIS ENVOLVENTE DE DATOS (DATA ENVELOPMENT ANALYSIS) # eficiencia T (x,y)^T for dato in range(dmu): # Un nuevo programa LP objeto prob nombre "hola" dea = LpProblem("dea", LpMinimize) # Definicion de la variable eficiencia eficiencia = LpVariable("eficiencia", 0) # Definicion de la variable lambda lambd = LpVariable.matrix("lambda", range(dmu), 0) # Definicion de funcion Objetivo dea += eficiencia for item in range(entr): dea += eficiencia*datIn0[dato,item] - lpDot(list(datIn0[:,item]),lambd) >= 0 for item in range(sal): dea += -datOut0[dato,item] + lpDot(list(datOut0[:,item]),lambd) >= 0 dea += lpSum(lambd) == 1 dea.solve() solu0[dato,mes] = value(dea.objective) # eficiencia T (x,y)^T+1 for dato in range(dmu): datIn1[dato,:] = datIn0[dato,:] datOut1[dato,:] = datOut0[dato,:] # Un nuevo programa LP objeto prob nombre "hola" dea = LpProblem("dea", LpMinimize) # Definicion de la variable eficiencia eficiencia = LpVariable("eficiencia", 0) # Definicion de la variable lambda lambd = LpVariable.matrix("lambda", range(dmu), 0) # Definicion de funcion Objetivo dea += eficiencia for item in range(entr): dea += eficiencia*datIn0[dato,item] - lpDot(list(datIn1[:,item]),lambd) >= 0 for item in range(sal): dea += -datOut0[dato,item] + lpDot(list(datOut1[:,item]),lambd) >= 0 dea += lpSum(lambd) == 1 dea.solve() solu1[dato,mes] = value(dea.objective) # eficiencia T+1 (x,y)^T for dato in range(dmu): datIn0[dato,:] = datIn1[dato,:] datOut0[dato,:] = datOut1[dato,:] # Un nuevo programa LP objeto prob nombre "hola" dea = LpProblem("dea", LpMinimize) # Definicion de la variable eficiencia eficiencia = LpVariable("eficiencia", 0) # Definicion de la variable lambda lambd = LpVariable.matrix("lambda", range(dmu), 0) # Definicion de funcion Objetivo dea += eficiencia for item in range(entr): dea += eficiencia*datIn1[dato,item] - lpDot(list(datIn0[:,item]),lambd) >= 0 for item in range(sal): dea += -datOut1[dato,item] + lpDot(list(datOut0[:,item]),lambd) >= 0 dea += lpSum(lambd) == 1 dea.solve() solu2[dato,mes] = value(dea.objective) print 'VER ESTO !!!!!', value(dea.objective) #print 'NORMAL' #print solu0 #print 'E (T) - (x,y)^T+1' #print solu1 Guarda(solu2,'resultados1.txt') # INDICE DE MALMQUIST malmquist = mat(zeros((dmu,meses))) EC = mat(zeros((dmu,meses))) PT = mat(zeros((dmu,meses))) for mes in range(meses-1): for d in range(dmu): EC[d,mes+1] = solu0[d,mes+1]/solu0[d,mes] PT[d,mes+1] = ((solu1[d,mes]/solu0[d,mes +1])*(solu0[d,mes]/solu2[d,mes]))**(1.0/2.0) malmquist[d,mes+1] = EC[d,mes+1] * PT[d,mes+1] ######################################################################################### Hopefully that serves to you, greetings, Carlos > _______________________________________________ > Help-glpk mailing list > [email protected] > http://lists.gnu.org/mailman/listinfo/help-glpk _______________________________________________ Help-glpk mailing list [email protected] http://lists.gnu.org/mailman/listinfo/help-glpk
