Constrained optimization in R My data is categorized by two different parameters (say F having n groups and S having m groups) and I want to get a relationship between the two. For example $F = ${$f_1 , f_2 , f_3$} = {$ 10,10,5$} and $S = ${$s_1, s_2 , s_3, s_4$} = {$8,8,8,4$}. (Read $f_1$ has 10 elements in it, $s_1$ has 8 elements in it.)
The problem is to get a relationship like $f_1 = p_{11}.s_1 + p_{12}.s_2 + p_{13}.s_3$ for all groups in F. The exact relationship is not possible (due to constraint 1 below), so we have to find the most approximate solution. The constraints are: 


*

*Sum of all the numbers in F (10+10+5) < sum of all the numbers in S
(8+8+8+4). Note that, this is a default property of data and we don't have to apply this.

*Total contribution of any $s_i$ can not be more than 100%. That is if $s_1$ contribute 40% in all 3 equations its total contribution will be 120% which should not happen. 

*Contribution can not be negative. 
I framed this as an optimization problem as
Given:
$F = (f_1, f_2, ... f_n)$
$S = (s_1, s_2, ... s_m)$
Each $f_i$ can be represented as: $f_i = p_{i1} s_1 + p_{i2}s_2 + ... + p_{im}s_m + \epsilon_i$
That is: $P = ${$p_{ij}$} is a $n * m$ matrix of row vectors $P_1, P_2, ..., P_n$ then $f_i = P_i.S+ \epsilon_i$ 
where . is the dot product (element by element multiplication).
I have to minimize $\sum_{i} \epsilon_i^2$, with following constraints:


*

*$p_{11}  + p_{12}+  ... + p_{1m}<=1$ That is some of each row of
matrix P is less than 1. 

*$p_{ij} >=0 $ for all i,j
I am trying to use constrOptim in R to solve this problem, however I am getting stuck with following:


*

*Framing the constraints

*How to create a code in which I just pass F and S and it gives me the matrix P and error vector as output.
Below is my code:
F = c(10,10,5)
S = c(8,8,8,4)
n = length(F)
m = length(S)
P_init = matrix(rep(0,m*n),nrow=n, ncol=m)

loss_fun <- function(P){

    T = S*P #proportion matrix * S
    F2 = rowSums(T) # Predicted values of F
    E = F - F2  # Error
    return(sum(E*E))
}
x = loss_fun(P_init)
z = constrOptim(P_init,loss_fun,NULL,ui=c(rep(-1,n)), ci=rep(-1,m))

Since the constraints are defined by: ui %*% theta - ci >= 0, I believe in my case ui = {-1,-1,-1}, theta = P_init ( a 3 x 4 matrix) and ci = {-1,-1,-1,-1}. However I get the following if I run the code. Error in ui %\*% theta : non-conformable arguments. Is theta = P_init is not true? Or there is some other error. As this is just one approach, I can also explore other approaches like optim or any other function. 
 A: I'm not sure if I'm right but help constrOptim {stats} say theta: numeric (vector) starting value (of length p): must be in the feasible region, and your theta is a matrix. 
I think here is near to your idea (Input is a vector and the function makes it a matrix).
[EDITED] NOT NEAR because of difference constraints. Please refer to OP's ANSWER.
F = c(10,10,5)
S = c(8,8,8,4)
n = length(F)
m = length(S)
P_init2 = rep(0, 12)

loss_fun2 <- function(a){
  P <- matrix(a, nrow=n, ncol=m)
  T = S*P
  F2 = rowSums(T)
  E = F - F2
  return(sum(E*E))
}

x = loss_fun2(rep(0, 12))
x  # 225
z = constrOptim(P_init2, loss_fun2, NULL, ui = rep(-1, n * m), ci = -1)

A: Using answer given by @crayfish and a detailed answer on how to put construct constraints here, I was able to come up with a solution.
#define F and S here
F = c(10,10,5)
S = c(8,8,9,8,4)

#loss_fun: to be minimized
loss_fun <- function(A){
    P = matrix(A, nrow = n,ncol = m, byrow=T)
    T = S*P #proportion matrix * S
    F2 = rowSums(T) # Predicted values of F
    E = F - F2  # Error
    return(sum(E*E))
}

n = length(F)
m = length(S)
#Initial solution (theta)
P_init = c(rep(0.1,n*m))

# Creating Constraints

vi = matrix(rep(0,n*m*n*(m+1)),ncol = n*m, byrow = T) 
for (i in 1:n){
    for (j in 1:(n*m)) {
        if (j <= m *  i & j > m * (i-1)) vi[i,j] = -1
    }
}
for (i in (n+1):(n*m+n)) {
    for (j in 1:(n*m)) {

        if ((i-n) == j)  vi[i,j] = 1
    }
}

myci = c(rep(-1,n),rep(0,n*m))

# check if initial value is in feasible region
vi %*% P_init - myci 

#run the optimization module
z = constrOptim(P_init,loss_fun,NULL,ui=vi, ci=myci)

#result
P_final = matrix(z$par,nrow=n,byrow=T)

Since creation of constraints may not be very easy to understand from the code, here is a visualization of constraints when n = 3 and m = 4
-1. p11 + -1. p12 + -1. p13 + -1. p14                                                                               >= -1
                                      -1. p21 + -1. p22 + -1. p23 + -1. p24                                         >= -1
                                                                          -1. p31 + -1. p32 + -1. p33 + -1. p34     >= -1
    p11                                                                                                             >= 0    
              p12                                                                                                   >= 0
                        p13                                                                                         >= 0
                                  p14                                                                               >= 0    
                                          p21                                                                       >= 0
                                                   p22                                                              >= 0
                                                              p23                                                   >= 0
                                                                       p24                                          >= 0
                                                                              p31                                   >= 0
                                                                                        p32                         >= 0
                                                                                                  p33               >= 0
                                                                                                            p34     >= 0

