Section 6.2

Page 251 #2

Nutritional Requirements-A rancher has determined that the minimum weekly nutritional requirements for an average-sized horse include 40 lb of protein, 20 lb of carbohydrates, and 45 lb of roughage. These are obtained from the following sources in varying amounts at the prices indicated:

items protein carbs roughage cost
hay 0.5 2.0 5.0 1.8
oats 1.0 4.0 2.0 3.5
feeding blocks 2.0 0.5 1.0 0.4
hpc 6.0 1.0 2.5 1.0
horse requirement 40.0 20.0 45.0 0.0

Formulate a mathematical model to determine how to meet the minimum nutritional requirements at minimum cost.

Solution

Variables

\(x_1=\mathrm{hay~per~bale}\\\)

\(x_2=\mathrm{oats~per~sack}\\\)

\(x_3=\mathrm{feeding~blocks~per~block}\\\)

\(x_4=\mathrm{high~protein~concentrate~per~sack}\)

Objective

\(\textrm{Minimize}~1.8 x_1+3.5x_2+0.4x_3+1.0x_4\)

Constraints

\(0.5x_1+x_2+2x_3+6x_4 \ge40\textrm{(protein)}\\\)

\(2x_1+4x_2+0.5x_3+x_4\ge20\textrm{(carbs)}\\\)

\(5x_1+2x_2+x_3+2.5x_4\ge45\textrm{(roughage)}\)

Page 264 #6

Maximize

\({10x+35y}\)

Subjected to

\(8x+6y\le48\textrm{(board-feet of lumber)}\\\)

\(4x+y\le20\textrm{(hours of capacity)}\\\)

\(y\ge5\textrm{(demand)}\\\)

\(x,y \ge0\textrm{(nonnegativity)}\)

obj.func <- function(x,y) 10*x + 35*y
s1 <- obj.func(0,5)
s2 <- obj.func(0,8)
s3 <- obj.func(9/4,5)
mydf <- data.frame(points=seq(1:3),x=c(0,0,9/4),y=c(5,8,5),obj_func = c(s1,s2,s3))
points x y obj_func
1 0.00 5 175.0
2 0.00 8 280.0
3 2.25 5 197.5
Solution

The objective function is maximized at point (0,8) with a value of 280.

Graphs
x <- seq(0, 6, by=0.1)

data <- data.frame(x = x, board_feet = (48-8*x)/6, hours = 20-4*x, demand = 5)
tidy <- gather(data, variable,  "y", -x) 

ggplotly(ggplot(tidy, aes(x=x, y=y, color=variable)) + geom_line())
myint <- data.frame(x=c(0,0,2.25),y=c(8,5,5))
ggplot(data.frame(x=c(-5,10)),aes(x)) + 
  stat_function(fun=function(x) -4/3*x+8, geom="line", aes(col='y=-4/3x+8')) +
  stat_function(fun=function(x) -4*x + 20, geom="line", aes(col='y=-4x+20')) +
  stat_function(fun=function(x)5, geom="line", aes(col='y=5')) + 
  geom_vline(xintercept=0, aes(col= 'x=0')) + 
  geom_hline(yintercept= 0, aes(col='y=0')) + 
  geom_point(data=myint, aes(x,y)) + 
  annotate('text', x = 0, y = 9.2, label="(0, 8)", size=3 ) +
  annotate('text', x = 0, y = 3.8, label="(0, 5)", size=3 ) + 
  annotate('text', x = 2.25, y = 3.8, label="(9/4, 5)", size=3 )

Page 295 #3

Use the curve-fitting criterion to minimize the sum of the absolute deviations for the following models and data set:

  1. \(y = ax\)
  2. \(y = ax^2\)
  3. \(y = ax^3\)

Model

\(y=ax\)

x <- c(7,14,21,28,35,42)
y <- c(8,41,133,250,280,297)
df1 <- data.frame(x=x, y=y)
a<- 0
b<- 50
r<- ((1+sqrt(5))/2)-1
cv1 = a + (1 - r)*(b - a)
cv2 = a + r*(b -a)
cost <- function(a, df){
  sum(abs(df$y - a*(df$x)))
}
crit_val_1<- cost(cv1,df1)
crit_val_2<- cost(cv2,df1)
d1 <- data.frame(a=a,b=b,cv1=cv1,cv2=cv2,crit_val_1=crit_val_1,crit_val_2=crit_val_2)
t <- 0.02
d <- 1
while(d > t){
  if(crit_val_1<=crit_val_2){
    b <- cv2
    cv1 <- a+(1-r)*(b-a)
    cv2 <- a+r*(b-a)
    crit_val_1 <- cost(cv1,df1)
    crit_val_2 <- cost(cv2,df1)
    d1 <- rbind(d1,c(a,b,cv1,cv2,crit_val_1,crit_val_2))
    d <- b-a
  } else{
    a <- cv1
    cv1 <- a+(1-r)*(b-a)
    cv2 <- a+r*(b-a)
    crit_val_1 <- cost(cv1,df1)
    crit_val_2 <- cost(cv2,df1)
    d1 <- rbind(d1,c(a,b,cv1,cv2,crit_val_1,crit_val_2))
    d <- b-a
  }
}
l<- length(d1$a)
d1$a[l] 
## [1] 7.066307
Calculate the Sum of the Absolute Deviations of the Model
d1$b[l]
## [1] 7.080309
cv1_final <- (d1$a[l]+ d1$b[l] ) /2
cv1_final
## [1] 7.073308
crit_val_1_optimal<- cost(cv1_final,df1)
crit_val_1_optimal
## [1] 199.5395

Model

\(y=ax^2\)

df2 <- data.frame(x=x, y=y, x2=x^2)
a<- 0
b<- 50
r<- ((1+sqrt(5))/2)-1
cv1 = a + (1 - r)*(b - a)
cv2 = a + r*(b -a)
cost2 <- function(a, df){
  sum(abs(df$y - a*(df$x2) ))
}
crit_val_1<- cost2(cv1,df2)
crit_val_2<- cost2(cv2,df2)
d2 <- data.frame(a=a,b=b,cv1=cv1,cv2=cv2,crit_val_1=crit_val_1,crit_val_2=crit_val_2)
t <- 0.02
d <- 1
while(d > t){
  if(crit_val_1<=crit_val_2){
    b <- cv2
    cv1 <- a+(1-r)*(b-a)
    cv2 <- a+r*(b-a)
    crit_val_1 <- cost2(cv1,df2)
    crit_val_2 <- cost2(cv2,df2)
    d2 <- rbind(d2,c(a,b,cv1,cv2,crit_val_1,crit_val_2))
    d <- b-a
  } else{
    a <- cv1
    cv1 <- a+(1-r)*(b-a)
    cv2 <- a+r*(b-a)
    crit_val_1 <- cost2(cv1,df2)
    crit_val_2 <- cost2(cv2,df2)
    d2 <- rbind(d2,c(a,b,cv1,cv2,crit_val_1,crit_val_2))
    d <- b-a
  }
}
l<- length(d2$a)
d2$a[l] 
## [1] 0.2232466
Calculate the Sum of the Absolute Deviations of the Model
d2$b[l]
## [1] 0.2372483
cv2_final <- (d2$a[l]+ d2$b[l] ) /2
cv2_final
## [1] 0.2302474
crit_val_2_optimal<- cost2(cv2_final,df2)
crit_val_2_optimal
## [1] 219.5671

Model

\(y=ax^3\)

df3 <- data.frame(x=x, y=y, x3=x^3)
a<- 0
b<- 50
r<- ((1+sqrt(5))/2)-1
cv1 = a + (1 - r)*(b - a)
cv2 = a + r*(b -a)
cost3 <- function(a, df){
  sum(abs(df$y - a*(df$x3)))
}
crit_val_1<- cost3(cv1,df3)
crit_val_2<- cost3(cv2,df3)
d3 <- data.frame(a=a,b=b,cv1=cv1,cv2=cv2,crit_val_1=crit_val_1,crit_val_2=crit_val_2)
t <- 0.02
d <- 1
while(d > t){
  if(crit_val_1<=crit_val_2){
    b <- cv2
    cv1 <- a+(1-r)*(b-a)
    cv2 <- a+r*(b-a)
    crit_val_1 <- cost3(cv1,df3)
    crit_val_2 <- cost3(cv2,df3)
    d3 <- rbind(d3,c(a,b,cv1,cv2,crit_val_1,crit_val_2))
    d <- b-a
  } else{
    a <- cv1
    cv1 <- a+(1-r)*(b-a)
    cv2 <- a+r*(b-a)
    crit_val_1 <- cost3(cv1,df3)
    crit_val_2 <- cost3(cv2,df3)
    d3 <- rbind(d3,c(a,b,cv1,cv2,crit_val_1,crit_val_2))
    d <- b-a
  }
}
l<- length(d3$a)
d3$a[l] 
## [1] 0
Calculate the Sum of the Absolute Deviations of the Model
d3$b[l]
## [1] 0.01400168
cv3_final <- (d3$a[l]+ d3$b[l] ) /2
cv3_final
## [1] 0.00700084
crit_val_3_optimal <- cost3(cv3_final,df3)
crit_val_3_optimal 
## [1] 433.7104
Solution
## [1] "The Optimal Value for the Model y=ax = 199.53946343722"
## [1] "The optimal value for the model y=ax^2 = 219.567063518652"
## [1] "The optimal value for the model y=ax^3 = 433.710399871649"