Farm Planning¶
Reference¶
SAS/OR example: http://go.documentation.sas.com/?docsetId=ormpex&docsetTarget=ormpex_ex8_toc.htm&docsetVersion=15.1&locale=en
SAS/OR code for example: http://support.sas.com/documentation/onlinedoc/or/ex_code/151/mpex08.html
Model¶
import sasoptpy as so
import pandas as pd
def test(cas_conn):
m = so.Model(name='farm_planning', session=cas_conn)
# Input Data
cow_data_raw = []
for age in range(12):
if age < 2:
row = {'age': age,
'init_num_cows': 10,
'acres_needed': 2/3.0,
'annual_loss': 0.05,
'bullock_yield': 0,
'heifer_yield': 0,
'milk_revenue': 0,
'grain_req': 0,
'sugar_beet_req': 0,
'labour_req': 10,
'other_costs': 50}
else:
row = {'age': age,
'init_num_cows': 10,
'acres_needed': 1,
'annual_loss': 0.02,
'bullock_yield': 1.1/2,
'heifer_yield': 1.1/2,
'milk_revenue': 370,
'grain_req': 0.6,
'sugar_beet_req': 0.7,
'labour_req': 42,
'other_costs': 100}
cow_data_raw.append(row)
cow_data = pd.DataFrame(cow_data_raw).set_index(['age'])
grain_data = pd.DataFrame([
['group1', 20, 1.1],
['group2', 30, 0.9],
['group3', 20, 0.8],
['group4', 10, 0.65]
], columns=['group', 'acres', 'yield']).set_index(['group'])
num_years = 5
num_acres = 200
bullock_revenue = 30
heifer_revenue = 40
dairy_cow_selling_age = 12
dairy_cow_selling_revenue = 120
max_num_cows = 130
sugar_beet_yield = 1.5
grain_cost = 90
grain_revenue = 75
grain_labour_req = 4
grain_other_costs = 15
sugar_beet_cost = 70
sugar_beet_revenue = 58
sugar_beet_labour_req = 14
sugar_beet_other_costs = 10
nominal_labour_cost = 4000
nominal_labour_hours = 5500
excess_labour_cost = 1.2
capital_outlay_unit = 200
num_loan_years = 10
annual_interest_rate = 0.15
max_decrease_ratio = 0.50
max_increase_ratio = 0.75
# Sets
AGES = cow_data.index.tolist()
init_num_cows = cow_data['init_num_cows']
acres_needed = cow_data['acres_needed']
annual_loss = cow_data['annual_loss']
bullock_yield = cow_data['bullock_yield']
heifer_yield = cow_data['heifer_yield']
milk_revenue = cow_data['milk_revenue']
grain_req = cow_data['grain_req']
sugar_beet_req = cow_data['sugar_beet_req']
cow_labour_req = cow_data['labour_req']
cow_other_costs = cow_data['other_costs']
YEARS = list(range(1, num_years+1))
YEARS0 = [0] + YEARS
# Variables
numCows = m.add_variables(AGES + [dairy_cow_selling_age], YEARS0, lb=0,
name='numCows')
for age in AGES:
numCows[age, 0].set_bounds(lb=init_num_cows[age],
ub=init_num_cows[age])
numCows[dairy_cow_selling_age, 0].set_bounds(lb=0, ub=0)
numBullocksSold = m.add_variables(YEARS, lb=0, name='numBullocksSold')
numHeifersSold = m.add_variables(YEARS, lb=0, name='numHeifersSold')
GROUPS = grain_data.index.tolist()
acres = grain_data['acres']
grain_yield = grain_data['yield']
grainAcres = m.add_variables(GROUPS, YEARS, lb=0, name='grainAcres')
for group in GROUPS:
for year in YEARS:
grainAcres[group, year].set_bounds(ub=acres[group])
grainBought = m.add_variables(YEARS, lb=0, name='grainBought')
grainSold = m.add_variables(YEARS, lb=0, name='grainSold')
sugarBeetAcres = m.add_variables(YEARS, lb=0, name='sugarBeetAcres')
sugarBeetBought = m.add_variables(YEARS, lb=0, name='sugarBeetBought')
sugarBeetSold = m.add_variables(YEARS, lb=0, name='sugarBeetSold')
numExcessLabourHours = m.add_variables(YEARS, lb=0,
name='numExcessLabourHours')
capitalOutlay = m.add_variables(YEARS, lb=0, name='capitalOutlay')
yearly_loan_payment = (annual_interest_rate * capital_outlay_unit) /\
(1 - (1+annual_interest_rate)**(-num_loan_years))
# Objective function
revenue = {year:
bullock_revenue * numBullocksSold[year] +
heifer_revenue * numHeifersSold[year] +
dairy_cow_selling_revenue * numCows[dairy_cow_selling_age,
year] +
so.expr_sum(milk_revenue[age] * numCows[age, year]
for age in AGES) +
grain_revenue * grainSold[year] +
sugar_beet_revenue * sugarBeetSold[year]
for year in YEARS}
cost = {year:
grain_cost * grainBought[year] +
sugar_beet_cost * sugarBeetBought[year] +
nominal_labour_cost +
excess_labour_cost * numExcessLabourHours[year] +
so.expr_sum(cow_other_costs[age] * numCows[age, year]
for age in AGES) +
so.expr_sum(grain_other_costs * grainAcres[group, year]
for group in GROUPS) +
sugar_beet_other_costs * sugarBeetAcres[year] +
so.expr_sum(yearly_loan_payment * capitalOutlay[y]
for y in YEARS if y <= year)
for year in YEARS}
profit = {year: revenue[year] - cost[year] for year in YEARS}
totalProfit = so.expr_sum(profit[year] -
yearly_loan_payment * (num_years - 1 + year) *
capitalOutlay[year] for year in YEARS)
m.set_objective(totalProfit, sense=so.MAX, name='totalProfit')
# Constraints
m.add_constraints((
so.expr_sum(acres_needed[age] * numCows[age, year] for age in AGES) +
so.expr_sum(grainAcres[group, year] for group in GROUPS) +
sugarBeetAcres[year] <= num_acres
for year in YEARS), name='num_acres')
m.add_constraints((
numCows[age+1, year+1] == (1-annual_loss[age]) * numCows[age, year]
for age in AGES if age != dairy_cow_selling_age
for year in YEARS0 if year != num_years), name='aging')
m.add_constraints((
numBullocksSold[year] == so.expr_sum(
bullock_yield[age] * numCows[age, year] for age in AGES)
for year in YEARS), name='numBullocksSold_def')
m.add_constraints((
numCows[0, year] == so.expr_sum(
heifer_yield[age] * numCows[age, year]
for age in AGES) - numHeifersSold[year]
for year in YEARS), name='numHeifersSold_def')
m.add_constraints((
so.expr_sum(numCows[age, year] for age in AGES) <= max_num_cows +
so.expr_sum(capitalOutlay[y] for y in YEARS if y <= year)
for year in YEARS), name='max_num_cows_def')
grainGrown = {(group, year): grain_yield[group] * grainAcres[group, year]
for group in GROUPS for year in YEARS}
m.add_constraints((
so.expr_sum(grain_req[age] * numCows[age, year] for age in AGES) <=
so.expr_sum(grainGrown[group, year] for group in GROUPS)
+ grainBought[year] - grainSold[year]
for year in YEARS), name='grain_req_def')
sugarBeetGrown = {(year): sugar_beet_yield * sugarBeetAcres[year]
for year in YEARS}
m.add_constraints((
so.expr_sum(sugar_beet_req[age] * numCows[age, year] for age in AGES)
<=
sugarBeetGrown[year] + sugarBeetBought[year] - sugarBeetSold[year]
for year in YEARS), name='sugar_beet_req_def')
m.add_constraints((
so.expr_sum(cow_labour_req[age] * numCows[age, year]
for age in AGES) +
so.expr_sum(grain_labour_req * grainAcres[group, year]
for group in GROUPS) +
sugar_beet_labour_req * sugarBeetAcres[year] <=
nominal_labour_hours + numExcessLabourHours[year]
for year in YEARS), name='labour_req_def')
m.add_constraints((profit[year] >= 0 for year in YEARS), name='cash_flow')
m.add_constraint(so.expr_sum(numCows[age, num_years] for age in AGES
if age >= 2) /
sum(init_num_cows[age] for age in AGES if age >= 2) ==
[1-max_decrease_ratio, 1+max_increase_ratio],
name='final_dairy_cows_range')
res = m.solve()
if res is not None:
so.pd.display_all()
print(so.get_solution_table(numCows))
revenue_df = so.dict_to_frame(revenue, cols=['revenue'])
cost_df = so.dict_to_frame(cost, cols=['cost'])
profit_df = so.dict_to_frame(profit, cols=['profit'])
print(so.get_solution_table(numBullocksSold, numHeifersSold,
capitalOutlay, numExcessLabourHours,
revenue_df, cost_df, profit_df))
gg_df = so.dict_to_frame(grainGrown, cols=['grainGrown'])
print(so.get_solution_table(grainAcres, gg_df))
sbg_df = so.dict_to_frame(sugarBeetGrown, cols=['sugerBeetGrown'])
print(so.get_solution_table(
grainBought, grainSold, sugarBeetAcres,
sbg_df, sugarBeetBought, sugarBeetSold))
num_acres = m.get_constraint('num_acres')
na_df = num_acres.get_expressions()
max_num_cows_con = m.get_constraint('max_num_cows_def')
mnc_df = max_num_cows_con.get_expressions()
print(so.get_solution_table(na_df, mnc_df))
return m.get_objective_value()
Output¶
In [1]: import os
In [2]: hostname = os.getenv('CASHOST')
In [3]: port = os.getenv('CASPORT')
In [4]: from swat import CAS
In [5]: cas_conn = CAS(hostname, port)
In [6]: import sasoptpy
In [7]: from examples.client_side.farm_planning import test
In [8]: test(cas_conn)
NOTE: Initialized model farm_planning.
NOTE: Added action set 'optimization'.
NOTE: Converting model farm_planning to OPTMODEL.
NOTE: Submitting OPTMODEL code to CAS server.
NOTE: Problem generation will use 8 threads.
NOTE: The problem has 143 variables (0 free, 13 fixed).
NOTE: The problem has 101 linear constraints (25 LE, 70 EQ, 5 GE, 1 range).
NOTE: The problem has 780 linear constraint coefficients.
NOTE: The problem has 0 nonlinear constraints (0 LE, 0 EQ, 0 GE, 0 range).
NOTE: The OPTMODEL presolver is disabled for linear problems.
NOTE: The LP presolver value AUTOMATIC is applied.
NOTE: The LP presolver time is 0.00 seconds.
NOTE: The LP presolver removed 86 variables and 71 constraints.
NOTE: The LP presolver removed 551 constraint coefficients.
NOTE: The presolved problem has 57 variables, 30 constraints, and 229 constraint coefficients.
NOTE: The LP solver is called.
NOTE: The Dual Simplex algorithm is used.
Objective
Phase Iteration Value Time
D 1 1 4.195000E+02 0
D 2 37 1.971960E+05 0
D 2 55 1.217192E+05 0
NOTE: Optimal.
NOTE: Objective = 121719.17286.
NOTE: The Dual Simplex solve time is 0.02 seconds.
NOTE: The output table 'SOLUTION' in caslib 'CASUSER(casuser)' has 143 rows and 6 columns.
NOTE: The output table 'DUAL' in caslib 'CASUSER(casuser)' has 101 rows and 4 columns.
NOTE: The CAS table 'solutionSummary' in caslib 'CASUSER(casuser)' has 13 rows and 4 columns.
NOTE: The CAS table 'problemSummary' in caslib 'CASUSER(casuser)' has 18 rows and 4 columns.
numCows
(0, 0) 10.000000
(0, 1) 22.800000
(0, 2) 11.584427
(0, 3) 0.000000
(0, 4) 0.000000
(0, 5) 0.000000
(1, 0) 10.000000
(1, 1) 9.500000
(1, 2) 21.660000
(1, 3) 11.005205
(1, 4) 0.000000
(1, 5) 0.000000
(2, 0) 10.000000
(2, 1) 9.500000
(2, 2) 9.025000
(2, 3) 20.577000
(2, 4) 10.454945
(2, 5) 0.000000
(3, 0) 10.000000
(3, 1) 9.800000
(3, 2) 9.310000
(3, 3) 8.844500
(3, 4) 20.165460
(3, 5) 10.245846
(4, 0) 10.000000
(4, 1) 9.800000
(4, 2) 9.604000
(4, 3) 9.123800
(4, 4) 8.667610
(4, 5) 19.762151
(5, 0) 10.000000
(5, 1) 9.800000
(5, 2) 9.604000
(5, 3) 9.411920
(5, 4) 8.941324
(5, 5) 8.494258
(6, 0) 10.000000
(6, 1) 9.800000
(6, 2) 9.604000
(6, 3) 9.411920
(6, 4) 9.223682
(6, 5) 8.762498
(7, 0) 10.000000
(7, 1) 9.800000
(7, 2) 9.604000
(7, 3) 9.411920
(7, 4) 9.223682
(7, 5) 9.039208
(8, 0) 10.000000
(8, 1) 9.800000
(8, 2) 9.604000
(8, 3) 9.411920
(8, 4) 9.223682
(8, 5) 9.039208
(9, 0) 10.000000
(9, 1) 9.800000
(9, 2) 9.604000
(9, 3) 9.411920
(9, 4) 9.223682
(9, 5) 9.039208
(10, 0) 10.000000
(10, 1) 9.800000
(10, 2) 9.604000
(10, 3) 9.411920
(10, 4) 9.223682
(10, 5) 9.039208
(11, 0) 10.000000
(11, 1) 9.800000
(11, 2) 9.604000
(11, 3) 9.411920
(11, 4) 9.223682
(11, 5) 9.039208
(12, 0) 0.000000
(12, 1) 9.800000
(12, 2) 9.604000
(12, 3) 9.411920
(12, 4) 9.223682
(12, 5) 9.039208
numBullocksSold numHeifersSold capitalOutlay numExcessLabourHours \
1 53.735000 30.935000 0.0 0.0
2 52.341850 40.757423 0.0 0.0
3 57.435807 57.435807 0.0 0.0
4 56.964286 56.964286 0.0 0.0
5 50.853436 50.853436 0.0 0.0
revenue cost profit
1 41494.530000 19588.466667 21906.063333
2 41153.336497 19264.639818 21888.696679
3 45212.490308 19396.435208 25816.055100
4 45860.056078 19034.285714 26825.770363
5 42716.941438 17434.354053 25282.587385
grainAcres grainGrown
(group1, 1) 20.000000 22.000000
(group1, 2) 20.000000 22.000000
(group1, 3) 20.000000 22.000000
(group1, 4) 20.000000 22.000000
(group1, 5) 20.000000 22.000000
(group2, 1) 0.000000 0.000000
(group2, 2) 0.000000 0.000000
(group2, 3) 3.134152 2.820737
(group2, 4) 0.000000 0.000000
(group2, 5) 0.000000 0.000000
(group3, 1) 0.000000 0.000000
(group3, 2) 0.000000 0.000000
(group3, 3) 0.000000 0.000000
(group3, 4) 0.000000 0.000000
(group3, 5) 0.000000 0.000000
(group4, 1) 0.000000 0.000000
(group4, 2) 0.000000 0.000000
(group4, 3) 0.000000 0.000000
(group4, 4) 0.000000 0.000000
(group4, 5) 0.000000 0.000000
grainBought grainSold sugarBeetAcres sugerBeetGrown sugarBeetBought \
1 36.620000 0.0 60.766667 91.150000 0.0
2 35.100200 0.0 62.670049 94.005073 0.0
3 37.836507 0.0 65.100304 97.650456 0.0
4 40.142857 0.0 76.428571 114.642857 0.0
5 33.476475 0.0 87.539208 131.308812 0.0
sugarBeetSold
1 22.760000
2 27.388173
3 24.550338
4 42.142857
5 66.586258
num_acres max_num_cows_def
1 200.0 130.000000
2 200.0 128.411427
3 200.0 115.433945
4 200.0 103.571429
5 200.0 92.460792
Out[8]: 121719.17286133638