[R] Fwd: Matrix Constraints in R Optim
Jeff Newmiller
jdnewmil at dcn.davis.ca.us
Sat Jun 18 00:57:36 CEST 2016
Always reply-to-all to keep the mailing list informed... I don't do private consulting online.
You look like you have plain text worked out, but we still cannot run your code.
Read the Posting Guide and [1]. You cannot attach Excel files here (and many R users can't or won't do that anyway) so you need to figure out dput. You will know you have a reproducible example when you don't read from any files and the code produces the error you are asking about when run in a newly opened R session. I think you should (!) only need the dput of my.data.var, definition of Error.func, and the call to optim if you do it right.
[1] http://stackoverflow.com/questions/5963269/how-to-make-a-great-r-reproducible-example
--
Sent from my phone. Please excuse my brevity.
On June 17, 2016 3:10:41 PM PDT, Narendra Modi <bjpmodi2016 at gmail.com> wrote:
>how to do that Jeff? I am newbie to R.
>
>I am posting the whole message again here and I made sure it is in
>plain-text format.
>
>
>
>
>
>
>my.data.matrix.inj <- as.matrix(my.data) #convert DATA FRAME to MATRIX
>my.data.matrix.inj
>
>
>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
>
>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
>
>
> # 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
>
>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),control=list(maxit
>=1000))
>sols
>
>
>
>
>
>
>On Fri, Jun 17, 2016 at 4:55 PM, Jeff Newmiller
><jdnewmil at dcn.davis.ca.us> wrote:
>>
>> Your code is corrupt because you failed to send your email in plain
>text format.
>>
>> You also don't appear to have all data needed to reproduce the
>problem. Use the dput function to generate R code form of a sample of
>your data.
>> --
>> Sent from my phone. Please excuse my brevity.
>>
>> On June 17, 2016 1:07:21 PM PDT, Priyank Dwivedi
><dpriyank23 at gmail.com> wrote:
>> >By mistake, I sent it earlier to the wrong address.
>> >
>> >---------- Forwarded message ----------
>> >From: Priyank Dwivedi <dpriyank23 at gmail.com>
>> >Date: 17 June 2016 at 14:50
>> >Subject: Matrix Constraints in R Optim
>> >To: r-help-owner at r-project.org
>> >
>> >
>> >Hi,
>> >
>> >Below is the code snippet I wrote in R:
>> >
>> >The basic idea is to minimize error by optimizing set of values (in
>> >this
>> >scenario 12) in the form of a matrix. I defined the matrix elements
>as
>> >vector "*my.data.var" * and then stacked it into a matrix called
>> >"*my.data.var.mat"
>> >in the error function. *
>> >
>> >The only part that I can't figure out is "what if the column sum in
>> >the *my.data.var.mat
>> >needs to be <=1"; that's the constraint/s.. Where do I introduce it
>in
>> >the
>> >OPTIM solver or elsewhere?*
>> >
>> >
>> >
>> >
>> >
>> >
>> >*my.data.matrix.inj* <- as.matrix(my.data) #convert DATA FRAME to
>> >MATRIX
>> >my.data.matrix.inj
>> >
>> >
>> >*my.data.matrix.time* <- as.matrix(my.data.2) #convert DATA FRAME
>to
>> >MATRIX
>> >my.data.matrix.time
>> >
>> >
>> >*my.data.matrix.prod* <- as.matrix(my.data) #convert DATA FRAME to
>> >MATRIX
>> >my.data.matrix.prod
>> >
>> >
>> >*my.data.var* <-
>>
>>c(2,0.8,0.5,0.2,0.2,0.1,10,0.01,0.02,0.2,0.1,0.01,2,0.8,0.5,0.2,0.2,0.1,10,0.01,0.02,0.2,0.1,0.01,2,0.8,0.5,0.2,0.2,0.1,10,0.01,0.02,0.2,0.1,0.01)
>> >my.data.var
>> >
>> >*my.data.qo* <- c(5990,150,199,996) #Pre-Waterflood Production
>> >
>> >*my.data.timet0* <- 0 # starting condition for time
>> >
>> >
>> >*#FUNCTIONQjk.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
>> >
>> >*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))
>> >
>> >
>> >
>> >--
>> >Best Regards,
>> >PD
>> >
>> >
>> >
>> >--
>> >Best Regards,
>> >Priyank Dwivedi
>> >
>> > [[alternative HTML version deleted]]
>> >
>> >______________________________________________
>> >R-help at r-project.org mailing list -- To UNSUBSCRIBE and more, see
>> >https://stat.ethz.ch/mailman/listinfo/r-help
>> >PLEASE do read the posting guide
>> >http://www.R-project.org/posting-guide.html
>> >and provide commented, minimal, self-contained, reproducible code.
>>
>> [[alternative HTML version deleted]]
>>
>> ______________________________________________
>> R-help at r-project.org mailing list -- To UNSUBSCRIBE and more, see
>> https://stat.ethz.ch/mailman/listinfo/r-help
>> PLEASE do read the posting guide
>http://www.R-project.org/posting-guide.html
>> and provide commented, minimal, self-contained, reproducible code.
[[alternative HTML version deleted]]
More information about the R-help
mailing list