% https://doi.org/10.1007/978-3-030-85165-1_18
% 4.4 n-site Phosphorylation

off rlabout;



rlset r;


{}


algebraic procedure varietyStarIsNotEmpty(xl, sfl, theo);
   rlsimpl(ex(xl, notzero xl and zero sfl), theo);


varietystarisnotempty


algebraic procedure varietyStarIsCosetInv(xl, sfl, theo);
   % all({g, x}, g<>0 and x<>0 and y<>0 and F(g)=0 and F(g*x)=0 impl F(g/x)=0);
   begin
      scalar il, gl, gxl, u1in;
      scalar x_in_variety, g_in_variety, gx_in_variety, g_over_x_in_variety;
      x_in_variety := zero sfl;
      il := for each x in xl collect getIndex(x, 'x);
      gl := for each i in il collect mkid(g, i);
      gxl := append(gl, xl);
      g_in_variety := sub(
         for each i in il collect mkid(x, i) = mkid(g, i),
         x_in_variety);
      gx_in_variety := sub(
         for each i in il collect mkid(x, i) = mkid(g, i) * mkid(x, i),
         x_in_variety);
      g_over_x_in_variety := sub(
         for each i in il collect mkid(x, i) = mkid(g, i) / mkid(x, i),
         x_in_variety);
      % Inverse:
      u1in := all(gxl, notzero gxl and g_in_variety and gx_in_variety
         impl g_over_x_in_variety);
      return u1in
   end;


varietystariscosetinv


algebraic procedure varietyStarIsCosetMult(xl, sfl, theo);
   % all({g, x, y}, g<>0 and x<>0 and y<>0 and F(g)=0 and F(g*x)=0 and F(g*y)=0 impl F(g*x*y)=0);
   begin
      scalar il, yl, gl, gxl, gxyl, u1, u2, g_in_variety, u1in, u2in;
      scalar x_in_variety, gx_in_variety, gy_in_variety, gxy_in_variety, g_over_x_in_variety;
      x_in_variety := zero sfl;
      il := for each x in xl collect getIndex(x, 'x);
      gl := for each i in il collect mkid(g, i);
      yl := for each i in il collect mkid(y, i);
      gxl := append(gl, xl);
      gxyl := append(gxl, yl);
      g_in_variety := sub(
         for each i in il collect mkid(x, i) = mkid(g, i),
         x_in_variety);
      gx_in_variety := sub(
         for each i in il collect mkid(x, i) = mkid(g, i) * mkid(x, i),
         x_in_variety);
      gy_in_variety := sub(
         for each i in il collect mkid(x, i) = mkid(g, i) * mkid(y, i),
         x_in_variety);
      gxy_in_variety := sub(
         for each i in il collect
            mkid(x, i) = mkid(g, i) * mkid(x, i) * mkid(y, i),
         x_in_variety);
      u2in := all(gxyl, notzero gxyl and g_in_variety and gx_in_variety and
         gy_in_variety impl gxy_in_variety);
      return u2in
   end;


varietystariscosetmult


symbolic operator getIndex;



symbolic procedure getIndex(v, base);
   begin scalar w;
      w := explode v;
      for i := 1:length explode base do
         w := cdr w;
      return compress w
   end;


getindex


algebraic procedure notzero(vl);
   rlsimpl for each v in vl  mkand v <> 0;


notzero


algebraic procedure zero(vl);
   rlsimpl for each v in vl  mkand v = 0;


zero


algebraic procedure rename(vl, fl);
   begin
      integer j;
      scalar sb, xj, xl;
      j := 0;
      xl := sb := {};
      for each v in vl do <<
         j := j + 1;
         xj := mkid(x, j);
         xl := cons(xj, xl);
         sb := cons(v = xj, sb)
      >>;
      return {reverse xl, sub(sb, fl), reverse sb}
   end;


rename


off allfac;


on rlnzden;


on rlqeaprecise;


off rlqedfs;


off rlqefbmma;


off rlqefbqepcad;


off rlqefbslfq;


on rlsiexpl;


on rlsiexpla;



procedure npho_vectorfield(n);
   begin scalar res, w;
      res := {-k_on0 * s0 * e + k_off0 * c0 + l_cat0 * d1};
      w := for i := 1:n-1 collect <<
         - mkid(k_on, i) * mkid(s, i) * e
         + mkid(k_off, i) * mkid(c, i)
         + mkid(k_cat, i-1) * mkid(c, i-1)
         - mkid(l_on, i-1) * mkid(s, i) * f
         + mkid(l_off, i-1) * mkid(d, i)
         + mkid(l_cat, i) * mkid(d, i+1)
      >>;
      res := append(res, w);
      w := for i := 0:n-1 collect <<
         mkid(k_on, i) * mkid(s, i) * e
         - (mkid(k_off, i) + mkid(k_cat, i)) * mkid(c, i)
      >>;
      res := append(res, w);
      w := for i := 1:n collect <<
         mkid(l_on, i-1) * mkid(s, i) * f
         - (mkid(l_off, i-1) + mkid(l_cat, i-1)) * mkid(d, i)
      >>;
      res := append(res, w);
      return res
   end;


npho_vectorfield


procedure npho_parameters(n);
   begin scalar res, w;
      res := for i := 0:n-1 collect
         mkid(k_on, i);
      w := for i := 0:n-1 collect
         mkid(k_off, i);
      res := append(res, w);
      w := for i := 0:n-1 collect
         mkid(k_cat, i);
      res := append(res, w);
      w := for i := 0:n-1 collect
         mkid(l_on, i);
      res := append(res, w);
      w := for i := 0:n-1 collect
         mkid(l_off, i);
      res := append(res, w);
      w := for i := 0:n-1 collect
         mkid(l_cat, i);
      res := append(res, w);
      return res
   end;


npho_parameters


procedure npho_variables(n);
   begin scalar res, w;
      res := for i := 0:n collect
         mkid(s, i);
      w := for i := 0:n-1 collect
         mkid(c, i);
      res := append(res, w);
      w := for i := 1:n collect
         mkid(d, i);
      res := append(res, w);
      res := append(res, {e, f});
      return res
   end;


npho_variables


n := 1;


n := 1


kl := npho_parameters n;


kl := {k_on0,

       k_off0,

       k_cat0,

       l_on0,

       l_off0,

       l_cat0}

vl := npho_variables n;


vl := {s0,

       s1,

       c0,

       d1,

       e,

       f}

ofl := npho_vectorfield n;


ofl := {c0*k_off0 + d1*l_cat0 - e*k_on0*s0,

         - c0*k_cat0 - c0*k_off0 + e*k_on0*s0,

         - d1*l_cat0 - d1*l_off0 + f*l_on0*s1}


w := rename(vl, ofl)$


xl := first w;


xl := {x1,x2,x3,x4,x5,x6}

fl := second w;


fl := {k_off0*x3 - k_on0*x1*x5 + l_cat0*x4,

        - k_cat0*x3 - k_off0*x3 + k_on0*x1*x5,

        - l_cat0*x4 - l_off0*x4 + l_on0*x2*x6}

third w;


{s0 = x1,s1 = x2,c0 = x3,d1 = x4,e = x5,f = x6}


pos := for each k in kl collect k > 0$



% Non-emptyness of V_k^*:
phi1 := varietyStarIsNotEmpty(xl, fl, pos);


phi1 := ex x1 ex x2 ex x3 ex x4 ex x5 ex x6 (x6 <> 0 and x5 <> 0 and x4 <> 0

 and x3 <> 0 and x2 <> 0 and x1 <> 0 and l_cat0*x4 + l_off0*x4 - l_on0*x2*x6 = 0

 and k_off0*x3 - k_on0*x1*x5 + l_cat0*x4 = 0

 and k_cat0*x3 + k_off0*x3 - k_on0*x1*x5 = 0)

w := rlqea(phi1, pos);


w := {{true,

       {x1 = -1,

        x2 = 1,

        x3 = 1,

              k_cat0
        x4 = --------,
              l_cat0

               - k_cat0 - k_off0
        x5 = --------------------,
                    k_on0

              k_cat0*l_cat0 + k_cat0*l_off0
        x6 = -------------------------------}}}
                      l_cat0*l_on0

phi1_ := rlsimpl(for each case in w mkor first case, pos);


phi1_ := true


% Shifted completeness under inverses:
phi2 := varietyStarIsCosetInv(xl, fl, pos);


phi2 := all g1 all g2 all g3 all g4 all g5 all g6 all x1 all x2 all x3 all x4 

all x5 all x6 (((x6 <> 0 and x5 <> 0 and x4 <> 0 and x3 <> 0 and x2 <> 0

 and x1 <> 0 and g6 <> 0 and g5 <> 0 and g4 <> 0 and g3 <> 0 and g2 <> 0

 and g1 <> 0) and ( - g2*g6*l_on0 + g4*l_cat0 + g4*l_off0 = 0

 and  - g1*g5*k_on0 + g3*k_off0 + g4*l_cat0 = 0

 and  - g1*g5*k_on0 + g3*k_cat0 + g3*k_off0 = 0) and (

 - g2*g6*l_on0*x2*x6 + g4*l_cat0*x4 + g4*l_off0*x4 = 0

 and  - g1*g5*k_on0*x1*x5 + g3*k_off0*x3 + g4*l_cat0*x4 = 0

 and  - g1*g5*k_on0*x1*x5 + g3*k_cat0*x3 + g3*k_off0*x3 = 0)) impl (

 - g2*g6*l_on0*x4 + g4*l_cat0*x2*x6 + g4*l_off0*x2*x6 = 0

 and  - g1*g5*k_on0*x3*x4 + g3*k_off0*x1*x4*x5 + g4*l_cat0*x1*x3*x5 = 0

 and  - g1*g5*k_on0*x3 + g3*k_cat0*x1*x5 + g3*k_off0*x1*x5 = 0))

phi2_ := rlqe(phi2, pos);


phi2_ := true


% Shifted completeness under multiplication:
phi3 := varietyStarIsCosetMult(xl, fl, pos);


phi3 := all g1 all g2 all g3 all g4 all g5 all g6 all x1 all x2 all x3 all x4 

all x5 all x6 all y1 all y2 all y3 all y4 all y5 all y6 (((y6 <> 0 and y5 <> 0

 and y4 <> 0 and y3 <> 0 and y2 <> 0 and y1 <> 0 and x6 <> 0 and x5 <> 0

 and x4 <> 0 and x3 <> 0 and x2 <> 0 and x1 <> 0 and g6 <> 0 and g5 <> 0

 and g4 <> 0 and g3 <> 0 and g2 <> 0 and g1 <> 0) and (

 - g2*g6*l_on0 + g4*l_cat0 + g4*l_off0 = 0

 and  - g1*g5*k_on0 + g3*k_off0 + g4*l_cat0 = 0

 and  - g1*g5*k_on0 + g3*k_cat0 + g3*k_off0 = 0) and (

 - g2*g6*l_on0*x2*x6 + g4*l_cat0*x4 + g4*l_off0*x4 = 0

 and  - g1*g5*k_on0*x1*x5 + g3*k_off0*x3 + g4*l_cat0*x4 = 0

 and  - g1*g5*k_on0*x1*x5 + g3*k_cat0*x3 + g3*k_off0*x3 = 0) and (

 - g2*g6*l_on0*y2*y6 + g4*l_cat0*y4 + g4*l_off0*y4 = 0

 and  - g1*g5*k_on0*y1*y5 + g3*k_off0*y3 + g4*l_cat0*y4 = 0

 and  - g1*g5*k_on0*y1*y5 + g3*k_cat0*y3 + g3*k_off0*y3 = 0)) impl (

 - g2*g6*l_on0*x2*x6*y2*y6 + g4*l_cat0*x4*y4 + g4*l_off0*x4*y4 = 0

 and  - g1*g5*k_on0*x1*x5*y1*y5 + g3*k_off0*x3*y3 + g4*l_cat0*x4*y4 = 0

 and  - g1*g5*k_on0*x1*x5*y1*y5 + g3*k_cat0*x3*y3 + g3*k_off0*x3*y3 = 0))

phi3_ := rlqe(phi3, pos);


phi3_ := true


coset := rlgsn(phi1_ and phi2_ and phi3_, form=dnf);


coset := true


phi4 := rlsimpl(sub(for each x in xl collect x=1, zero fl), pos);


phi4 := l_cat0 + l_off0 - l_on0 = 0 and k_off0 - k_on0 + l_cat0 = 0

 and k_cat0 + k_off0 - k_on0 = 0


group := rlgsn(coset and phi4, form=dnf);


group := k_off0 - k_on0 - l_off0 + l_on0 = 0 and k_off0 - k_on0 + l_cat0 = 0

 and k_cat0 + k_off0 - k_on0 = 0


end;

