Skip to content
Prev 361856 / 398506 Next

Matrix Constraints in R Optim

All,
Here are the dput files of the input data to the code.

Thanks for any advice.

I am adding the entire code below too just in case.


file <- file.path("Learning R","CRM_R_Ver4.xlsx")
file
my.data <- readWorksheetFromFile(file,sheet=1,startRow=1)
str(my.data)  # DATA FRAME
my.data.matrix.inj <- as.matrix(my.data)  #convert DATA FRAME to MATRIX
my.data.matrix.inj

dput(my.data.matrix.inj,"my.data.matrix.inj.txt")


my.data.2 <- readWorksheetFromFile(file,sheet=2,startRow=1)
str(my.data.2)  # DATA FRAME
my.data.matrix.time <- as.matrix(my.data.2)  #convert DATA FRAME to MATRIX
my.data.matrix.time

dput(my.data.matrix.time,"my.data.matrix.time.txt")

my.data <- readWorksheetFromFile(file,sheet=3,startRow=1)
str(my.data)  # DATA FRAME
my.data.matrix.prod <- as.matrix(my.data)  #convert DATA FRAME to MATRIX
my.data.matrix.prod

dput(my.data.matrix.prod,"my.data.matrix.prod.txt")

 # my.data.var <- vector("numeric",length = 24)
 # my.data.var

my.data.var <- c(10,0.25,0.25,0.25,0.25,0.25,
                 10,0.25,0.25,0.25,0.25,0.25,
                 10,0.25,0.25,0.25,0.25,0.25,
                 10,0.25,0.25,0.25,0.25,0.25)
my.data.var

dput(my.data.var,"my.data.var.txt")


my.data.qo <- c(5990,150,199,996)   #Pre-Waterflood Production
my.data.timet0 <- 0 # starting condition for time

#FUNCTION
Qjk.Cal.func <- function(my.data.timet0,my.data.qo,my.data.matrix.time,
                         my.data.matrix.inj,
my.data.matrix.prod,my.data.var,my.data.var.mat)
{

  qjk.cal.matrix <- matrix(,nrow = nrow(my.data.matrix.prod),
ncol=ncol(my.data.matrix.prod))

  count <- 1
  number <- 1
  for(colnum in 1:ncol(my.data.matrix.prod))   # loop through all PROD
wells columns
  {
    sum <-0
    for(row in 1:nrow(my.data.matrix.prod)) #loop through all the rows
    {
      sum <-0
      deltaT <-0
      expo <-0


        for(column in 1:ncol(my.data.matrix.inj)) #loop through all
the injector columns to get the PRODUCT SUM
         {
            sum = sum +
my.data.matrix.inj[row,column]*my.data.var.mat[colnum,number+column]
         }

      if(count<2)
      {
        deltaT<- my.data.matrix.time[row]
      }
      else
      {deltaT <- my.data.matrix.time[row]-my.data.matrix.time[row-1]}


      expo <- exp(-deltaT/my.data.var.mat[colnum,1])
# change here too

      if(count<2)
      {
        qjk.cal.matrix[row,colnum] = my.data.qo[colnum]*expo + (1-expo)*sum
      }
      else
      {
        qjk.cal.matrix[row,colnum]=qjk.cal.matrix[row-1,colnum]*expo +
(1-expo)*sum
      }
      count <- count+1
    }

    count <-1
  }

  qjk.cal.matrix      # RETURN CALCULATED MATRIX TO THE ERROR FUNCTION

}


# ERROR FUNCTION - FINDS DIFFERENCE BETWEEN CAL. MATRIX AND ORIGINAL
MATRIX. Miminize the Error by changing my.data.var

Error.func <- function(my.data.var)
{
  #First convert vector(my.data.var) to MATRIX aand send it to
calculate new MATRIX
  my.data.var.mat <- matrix(my.data.var,nrow =
ncol(my.data.matrix.prod),ncol = ncol(my.data.matrix.inj)+1,byrow =
TRUE)

  Calc.Qjk.Value <- Qjk.Cal.func(my.data.timet0,my.data.qo,my.data.matrix.time,
                                 my.data.matrix.inj,
my.data.matrix.prod,my.data.var,my.data.var.mat)


  diff.values <- my.data.matrix.prod-Calc.Qjk.Value    #FIND
DIFFERENCE BETWEEN CAL. MATRIX AND ORIGINAL MATRIX


  Error <- ((colSums ((diff.values^2), na.rm = FALSE, dims =
1))/nrow(my.data.matrix.inj))^0.5    #sum of square root of the diff
  print(paste(Error))

  Error_total <- sum(Error,na.rm=FALSE)/ncol(my.data.matrix.prod)   #
total avg error


  Error_total
}

# OPTIMIZE

sols<-optim(my.data.var,Error.func,method="L-BFGS-B",upper=c(Inf,1,1,1,1,1,Inf,1,1,1,1,1,Inf,1,1,1,1,1,Inf,1,1,1,1,1),
      lower=c(0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0))

sols
On 17 June 2016 at 16:55, Jeff Newmiller <jdnewmil at dcn.davis.ca.us> wrote: