Question: I am really struggling with this subs thing!

Working.mw

This one is wokring fine.

However,

See the new worksheet: testP.mw

with the source code: tmP.c

######################################################################
intM:=proc(K::integer, C::integer)

description "Creat interaction terms (matrix form).";

local
    M:=Matrix(K,C,symbol=mix),
    tolt,cols,rows;

# Set constraints

    M(..,1):=0;
    M(1..2,..):=0;
    
    tolt:=convert(M,`+`);
    cols:=MTM:-sum(M,1);
    rows:=MTM:-sum(M,2);
    
    M(2,..):=-cols;
    M(..,1):=-rows;
    M(2,1):=tolt;
    
    M;
    
end proc;



myM:=proc(
    K::integer,
    C::integer,
    tvarp::truefalse:=false,
    phiMix::truefalse:=false
    )

local
    p,     # Vectors of phi and p
    int0;  # Matrix of interaction terms

# take care of p

    if tvarp then
        if phiMix then
            int0:=intM(K,C);
            p:=Array(1..K,1..C,(j,c)->1/(1+exp(-(mu+tau[j]+eta[c]+int0[j,c]))));
            # p:=Array(1..K,1..C,(j,c)->mu+tau[j]+eta[c]+int0[j,c]);
        else
            if C=1 then
                p:=Array(1..K,1..C,(j,c)->1/(1+exp(-(mu+tau[j]))));
                # p:=Array(1..K,1..C,(j,c)->mu+tau[j]);
            else
                p:=Array(1..K,1..C,(j,c)->1/(1+exp(-(mu+tau[j]+eta[c]))));
                # p:=Array(1..K,1..C,(j,c)->mu+tau[j]+eta[c]);
            end if;
        end if;
    else
            p:=Array(1..K,1..C,(j,c)->1/(1+exp(-(mu+eta[c]))));
            # p:=Array(1..K,1..C,(j,c)->mu+eta[c]);
    end if;
    
    p(1,..):=0;
    p:=subs({tau[2]=0,eta[1]=0},p);
    
    p;
    
end proc;

myrand:=evalf[2](rand(10..90)/100):
######################################################################

I am really struggling why it does not replace the ans into the original matrix!

Many thanks,

casper