Creating an equation library (Updated 03-FEB-2017)

[Han]

(This post has been edited; the original content of the first post is at the bottom of this post)

MARCH 2017: This project now exists as an actual app. Please see http://hpmuseum.org/forum/thread-7725.html

------------

Version 0.063 now makes use of the SVD package I wrote and posted here: http://www.hpmuseum.org/forum/thread-497...67743.html

This version falls back on the SVD package I wrote for certain classes of matrices -- mostly those with singular values of 0. The Prime's implementation computes the left orthonormal vectors based on the reciprocal of the singular values. So when a singular value is zero, this results in undefined orthonormal vectors. There are still many errors that need to be trapped, but if you stall the SVD program in addition to the code below, you should be able to solve any system, now -- provided your guesses does not result in a point \( \mathbf{x} \) such that \( \mathbf{F}(\mathbf{x} )\) or its Jacobian is undefined.

I will eventually write some documentation, but for now the solve is a simple Newton solver. Briefly,

\[ \mathbf{F}(\mathbf{x}) =
\left[ \begin{array}{c}
f_1 (x_1, x_2, \dotsm, x_n)\\
f_2 (x_1, x_2, \dotsm, x_n) \\
\vdots \\
f_m (x_1, x_2, \dotsm, x_n)
\end{array} \right] \]

where the \( f_i \)'s are the equations set to 0. For example, \( x^2 + xy = 9 \) becomes \( f(x,y) = x^2+xy - 9 \). The linearization of \( \mathbf{F} \)
is
\[ \mathbf{F} (\mathbf{x}_{n+1} ) = \mathbf{F} (\mathbf{x}_n ) + J_\mathbf{F} \cdot \underbrace{(\mathbf{x}_{n+1} - \mathbf{x}_{n} ) }_{\Delta \mathbf{x}}\]
where \( J_{\mathbf{F}} \) is the Jacobian (evaluated at \( \mathbf{x}_{n} \) ). The solver solves for \( \Delta \mathbf{x} \) by computing the pseudo inverse of \( J_\mathbf{F} \) -- which is exactly the inverse if the Jacobian if it is non-singular. I won't get too much SVD and pseudo-inverse until I write up the documentation.

Code:

ssVersion:="Equation Library 0.063 by Han Duong";

// colors
ssFG:=#0h;
ssBG:=#FFFFFFh;
ssR:=#FF0000h;
ssG:=#FF00h;
ssB:=#FFh;


// error messages
ssInvDat:="Invalid library data.";
ssNoDat:="Unable to locate library data.";
ssNullName:="Null name not allowed.";
ssInvName:="Invalid name.";
ssAbort:="Operation aborted; improper initialization.";
ssBadVar:="Conflict with existing variable:\n";
ssSubmitBug:="This should never happen!\nPLease submit bug report.";
ssNotEqn:="Non-equation on page ";
ssNonVar:="No such var. on page ";
ssNoLibDat:="No library data found.\nCreating empty library and empty\ninitial system.";
ssEOW:="Overwrite disabled and name already exists:\n";
ssNotSaved:="SAVE CANCELED!";
ssMustChoose:="Select a system or create a new one.";
ssNoVars:="No variables for this system!";
ssNoSystems:="No systems in library!";
ssBadJac:="Jacobian not invertible; singular or non-square?";
ssFUndef:="Error evaluating F!";

// app messages
ssTNewVar:="Create/Edit Variable";
ssLName:="Name:";
ssLDescr:="Description:";
ssLOWrite:="Overwrite";
ssHName:="Enter the variable name";
ssHDescr:="Enter a description for the variable";
ssHOWrite:="Overwrite if variable already exists?";
ssTSaveSys:="Save Current System?";
ssHSysName:="Enter a name for the system.";
ssLEqns:="Eqns:";
ssLVars:="Vars:";
ssHEqns:="Verify equations in system.";
ssHVars:="Verify variables for system.";
ssTNewSys:="New System";
ssTChooseSys:="Choose a System";
ssHMakeConst:="Make this variable constant?";
ssTDeleteSys:="Select Systems to Delete";
ssHMarkDel:="Mark for deletion?";


ssTSettings:="Equation Library Settings";
ssLSaving:="When switching systems, always";
ssCSaving:={ "prompt to save", "save automatically", "discard all changes"};
ssHSaving:="Select default save behavior";


// misc
ssVLoadEqns:="Populating equations in Symb view...";
ssVLoadVars:="Setting up equation variables...";
ssVCheckVars:="Verifying variables in system...";
ssVCheckLibDat:="Verifying library data...";
ssVUpdateEqns:="Equations modified; updating lib data...";
ssVUpdateVars:="Variables modified; updating lib data...";
ssVNoEqns:="(No equations)";
ssVNoVars:="(No variables)";
ssVVarMade:="Created variable: ";


// View menu options
ssTSEP:="Select Equations Page";
ssTSelPage:="Select New Page";
ssTNewPage:="Create New Page";


ssMsgRow;
export ssInit;
ssLibFN;
export ssCurSys;
export ssCurSysIndex;
export ssCurLib;
ssLibSize;
export ssSysTitles;
ssVarKey:="SSVersion";
export ssEqPages;
export ssEqPage;
export ssNullSys;
export ssWarning;
export ssSolveVars:={};
export ssAutoVar:=1;


// DO NOT MODIFY BELOW HERE

ssGlobal:={ "A", "B", "C", "D", "E", "F", "G",
            "H", "I", "J", "K", "L", "M", "N",
            "O", "P", "R", "S", "T", "U", "V",
            "W", "X", "Y", "Z", "θ",
            "Z0", "Z1", "Z2", "Z3", "Z4",
            "Z5", "Z6", "Z7", "Z8", "Z9"
};


ssUpdateSys()
begin
  local n;

  ssEqPages:={1,10};
  n:=size(ssCurSys(2));
  ssNullSys:=NOT(n);
  ssEqPages(1):=ip((n-1)/10);
  ssEqPages(2):=(n-1) mod 10;
  ssEqPage:=0;
  for n from 1 to size(ssCurSys(4)) do
    if NOT(ssCurSys(7,n)) then ssSolveVars(0):=ssCurSys(4,n); end;
  end;
  DelAVars(AVars);
end;


ssInitSSApp()
begin
  local j,run:=1;
  local cmd:="";
  if ssInit then return(0); end;

  dimgrob_p(G1,320,240);

  ssLibFN:="Equation Library.lib";
  ssVarKey:="SSVersion";
  ssWarning:=1;

  iferr
    ssCurLib:=AFiles(ssLibFN);
  then
    msgbox(ssNoLibDat);
    ssCurLib:={};
    AFiles(ssLibFN):=ssCurLib;
  end;
  for j from 0 to 9 do
    cmd:="E" + string(j,1,0) + ":=" + string("");
    expr(cmd);
  end;
  // always start app with new system
  ssCurSys:={ssTNewSys, {}, {}, {}, {}, {}, {} };
  ssLibSize:=size(ssCurLib);
  if ssLibSize then
    ssSysTitles:=makelist(ssCurLib(j,1),j,1,ssLibSize);
  end;
  ssCurSysIndex:=0;
  ssUpdateSys();
  ssInit:=1;
end;


ssPrint(msg)
begin
  textout_p(msg,G1,0,12*ssMsgRow,1,ssFG,320,ssBG);
  line_p(G1,0,12*ssMsgRow+12,12,12*ssMsgRow+12,ssR);
  ssMsgRow:=(ssMsgRow+1) mod 20;
end;


ssError(msg)
begin
  textout_p(msg,G1,0,12*ssMsgRow,1,ssR,320,ssBG);
  line_p(G1,0,12*ssMsgRow+12,12,12*ssMsgRow+12,ssR);
  ssMsgRow:=(ssMsgRow+1) mod 20;
end;

ssWarn(msg)
begin
  if ssWarning then
    msgbox(msg);
  end;
end;


ssCheckType(list,t)
begin
  local n:=size(list);
  if n then
    if type(list(n))<>t then return(1); end;
  else
    return(1);
  end;
  return(0);
end;



//******************************************
// Library Data format
// { sys1, sys2, ... , sysN }
//
// Each sys (system) is a list of the form
// {
//   "Title",
//   { "eq1", "eq2", ... },
//   { sel1, sel2, ... },
//   { "var1", "var2", ... },
//   { "des1", "des2", ... },
//   { val1, val2, ... },
//   { con1, con2, ... }
// }
//
// eq  = equation
// sel = 0 or 1 for selection of corresponding
//       equation in the system
// var = variable
// des = description of variable
// val = value of corresponding variable
// con = 0 or 1 for whether var is constant
//
// returns 0 if error or 1 if OK
//******************************************
export ssCheckLibDat()
begin

  ssPrint(ssVCheckLibDat);
  ssInitSSApp();

  local syseqs,j,k;
  local lib, isbad:=0;
  local types:={2,2,0,2,2,0,0};

  iferr
    if ssLibSize then
      for j from 1 to ssLibSize do
        lib:=ssCurLib(j);
        if type(lib)<>6 then
          ssError(ssInvDat); return(0);
        else
          k:=size(lib);
          if (k<1) OR (k>7) then
            ssError(ssInvDat); return(0);
          end;
          for k from 2 to 7 do
            if type(lib(k))<>6 then
              ssError(ssInvDat); return(0);
            end;
            if ssCheckType(lib(k),types(k)) then
              ssError(ssInvDat); return(0);
            end;
          end; // for k
        end; // if type(lib)
      end; // for j
    end; // if ssLibSize
  then
    ssError(ssNoDat); return(0);
  end;
  return(1);
end;


// initialize variables
ssInitVars()
begin

  local varslist:=ssCurSys(4);
  local valslist:=ssCurSys(6);
  local n:=size(varslist);
  local j,varval;


  iferr
    if (AVars(ssVarKey) == ssVersion) then   
      return(0);
    end;
  then
    // nop
  end;

  ssPrint(ssVLoadVars);
  for j from 1 to n do
    if pos(ssGlobal,varslist(j)) then
      expr(varslist(j) + ":=" + valslist(j));
    else
      varval:=CAS(EVAL(varslist(j)));
      purge(EVAL(varslist(j)));
      if (type(CAS(EVAL(varslist(j)))) <> 8) then
        ssWarn(ssBadVar+varslist(j));
      end;
      if ( (type(varval) == 0) OR (type(varval) == 3) ) then
        AVars(varslist(j)):=varval;
      else
        AVars(varslist(j)):=valslist(j);
      end;
    end;
  end;
  AVars(ssVarKey):=ssVersion;

end;


// checks to see if AVars and system of equations
// vars match; excludes version key
export ssCheckVars()
begin

  ssPrint(ssVCheckVars);
  ssInitSSApp();
  ssInitVars();

  local j,n;
  local avarlist:=AVars();
  local varlist:=ssCurSys(4);
  n:=size(varlist);

  for j from 1 to n do
    if NOT(pos(avarlist,varlist(j)) OR pos(ssGlobal,varlist(j))) then
      return(0);
    end;
  end;

  n:=size(avarlist);
  for j from 1 to n do
    if ( (avarlist(j) <> ssVarKey) AND NOT(pos(varlist,avarlist(j))) ) then
      return(0);
    end;
  end;
  return(1);

end;


// do we have any equations
//  1: yes, equations in symb view
//  0: no equations at all
// -1: yes, but only saved in library
ssHaveEqns()
begin
  local j;
  for j from 0 to 9 do
    if ISCHECK(j) then return(1); end;
  end;

  if (ssLibSize == 0) then return(0); end;

  if (type(ssCurSys) <> 6) then return(0); end;

  if size(ssCurSys) then
    if size(ssCurSys(2)) then
      return(-1);
    else
      return(0);
    end;
  else
    return(0);
  end;
end;


// UI for creating new variables
ssCreateVar()
begin
  local run:=1;
  local n;
  local varname:="", varinfo:="";
  local ow:=0;

  ssInitSSApp();
  ssInitVars();

  while run do
    if input(
      {
        {varname, [2], {30, 65, 0}},
        {varinfo, [2], {30, 65, 1}},
        {ow, 1, {30,10,2}}
      },
      ssTNewVar,
      { ssLName, ssLDescr, ssLOWrite },
      { ssHName, ssHDescr, ssHOWrite }
    )
    then

      if size(varname)<1 then
        msgbox(ssNullName);
      else
        n:=pos(AVars,varname);
        if n then
          if ow then
            n:=pos(ssCurSys(4),varname);
            ssCurSys(5,n):=varinfo;
            run:=0;
          else
            msgbox(ssEOW + varname);
          end;
        else
          iferr

            if NOT(pos(ssGlobal,varname)) then
              AVars(varname):=0;
            end;

          then
            msgbox(ssInvName);
          else
            purge(EVAL(varname));
            ssCurSys(4,0):=varname;
            ssCurSys(5,0):=varinfo;
            ssCurSys(6,0):=0;
            ssCurSys(7,0):=0;
            msgbox(ssVVarMade + varname);
            varname:=""; varinfo:="";
          end;
        end;
      end; // end null name

    else
      run:=0;
    end;

  end; // end while
end;


// save current equations
ssSaveEqns()
begin
  local j,k,n1,n2;
  local cmd:="";
  local oldeqn:="";
  local neweqn:="";
  local markadd:=0;
  local neweqns:={};
  local teqnlist:={};
  local tsellist:={};

  if ssCheckVars() then

    ssPrint(ssVLoadEqns);
    n1:=ssEqPage*10+1;
    if (ssEqPage < ssEqPages(1)) then
      n2:=n1+9;
    else
      n2:=n1+ssEqPages(2);
    end;

    for j from n1 to n2 do
      k:= j mod 10;
      cmd:="string(E" + string(k,1,0) + ",1,12)";

      iferr
        neweqn:=expr(cmd);
      then
        neweqn:="";
      else
        if ssNullSys then
          oldeqn:="";
        else
          cmd:="E" + string(k,1,0) + ":=" + string(ssCurSys(2,j));
          iferr expr(cmd); then end;
          cmd:="string(E" + string(k,1,0) + ",1,12)";
          iferr oldeqn:=expr(cmd); then end;
        end;
        if (oldeqn <> neweqn) then
          if ssNullSys then
            ssCurSys(2,0):=neweqn;
          else
            ssCurSys(2,j):=neweqn;
          end;
        end;
      end;

      cmd:="E" + string(k,1,0) + ":=" + string(neweqn);
      iferr expr(cmd); then
        // msgbox(ssSubmitBug);
      end;

      if ISCHECK(k) then
        ssCurSys(3,j):=1;
      else
        ssCurSys(3,j):=0;
      end;

      if (size(neweqn) == 0) then
        ssCurSys(3,j):=-1; // mark for delete
      end;
    end;  // end for j

    n2:=n1+9;
    for k from j to n2 do
      n1:=k mod 10;
      cmd:="string(E" + string(n1,1,0) + ",1,12)";
      iferr
        neweqn:=expr(cmd);
      then
        // nop
      else
        if size(neweqn) then
          neweqns(0):=neweqn;
          neweqns(0):=ISCHECK(n1);
          markadd:=1;
        end;
      end;
    end; // end for k

    // rebuilt existing equations
    k:=0;
    n2:=size(ssCurSys(2));
    for j from 1 to n2 do
      if (ssCurSys(3,j) >= 0) then
        teqnlist(0):=ssCurSys(2,j);
        tsellist(0):=ssCurSys(3,j);
        k:=k+1;
      end;
    end;
    ssCurSys(2):=teqnlist;
    ssCurSys(3):=tsellist;

    if markadd then
      n2:=size(neweqns)/2;
      for j from 1 to n2 do
        ssCurSys(2,0):=neweqns(2*j-1);
        ssCurSys(3,0):=neweqns(2*j);
        k:=k+1;
      end;
    end;

    ssEqPages(1):=ip((k-1)/10);
    ssEqPages(2):=(k-1) mod 10;
    ssEqPage:=min(max(0,ssEqPage),ssEqPages(1));

    if k then
      ssNullSys:=0;
    else
      ssNullSys:=1;
    end;

    return(1);

  else
    msgbox(ssAbort);
    return(0);
  end;
end;


// set up equations in Symb view
export ssLoadSymb()
begin

  local j,k,n1,n2;
  local cmd:="";

  // no equations to populate Symb view?
  if ssNullSys then 
    for j from 0 to 9 do
      cmd:="E" + string(j,1,0) + ":=" + string("");
      expr(cmd);
    end;
    return(1);
  end;

  ssPrint(ssVLoadEqns);
  n1:=ssEqPage*10+1;
  if (ssEqPage < ssEqPages(1)) then
    n2:=n1+9;
  else
    n2:=n1+ssEqPages(2);
  end;
  for j from n1 to n2 do
    k:= j mod 10;
    cmd:="E" + string(k,1,0) + ":=" + string(ssCurSys(2,j));
    iferr expr(cmd); then end;
    if (ssCurSys(3,j)==1) then CHECK(k); else UNCHECK(k); end;
  end;

  n2:=n1+9;
  for k from j to n2 do
    n1:=k mod 10;
    cmd:="E" + string(n1,1,0) + ":=" + string("");
    expr(cmd);
  end;

  return(1);

end;


// select an equations page
ssSelectEqPage()
begin
  ssInitSSApp();
  local j;
  local n:=ssEqPages(1)+1;
  local page:=ssEqPage+1;
  local curpage:=page;
  local pagelist:={};
  local title:="[" + string(page,1,0) + "/" + string(ssEqPages(1)+1,1,0) + "] " + ssTSelPage;
  
  pagelist:=makelist(j,j,1,n);
  pagelist(0):=ssTNewPage;

  if choose(page,title,pagelist) then

    if (curpage <> page) then
      ssSaveEqns(); 

      if (page == (n+1)) then

        if (ssNullSys == 0) then 
          ssEqPages(1):=n+1;

          // fill last page
          if (ssEqPages(2) < 9) then
            n:=9-ssEqPages(2);
            for j from 1 to n do
              ssCurSys(2,0):="";
              ssCurSys(3,0):=-1; // mark for deletion
            end;
          end;
        end;

        // create new page
        for j from 1 to 10 do
          ssCurSys(2,0):="";
          ssCurSys(3,0):=-1; // mark for deletion
        end;

        n:=size(ssCurSys(2));
        ssEqPage:=ip((n-1)/10);
        ssEqPages(1):=ssEqPage;
        ssEqPages(2):=(n-1) mod 10;
        ssNullSys:=0;

      else

        // not creating new blank page
        ssEqPage:=min(max(0,page-1),ssEqPages(1));

      end; // if newpage or existing page

    end; // if selpage <> page

    ssLoadSymb();

  end; // if choose
end;


// set new system
export ssSetSystem(n)
begin
  if (n > ssLibSize) then
    return(0);
  end;

  ssCurSysIndex:=n;
  ssCurSys:=ssCurLib(n);
  ssUpdateSys();
  ssInitVars();
  ssLoadSymb();
end;


ssErrorPage(n,msg,eqn)
begin
  local p:=ip((n-1)/10)+1;
  msgbox(msg + n + ":\n" + eqn);
end;


// automatically add existing variable to vars list
// should only happen with global variables
ssAutoAddVar(var)
begin
  local val;
  local t;

  iferr
    val:=expr(var);
  then
    val:=0;
  end;

  t:=type(val);
  if ( (t<>0) OR (t<>3) ) then
    val:=0;
  end;   

  // convert existing vars to app var
  // exceptions: global system vars cannot
  // be app vars in current firmware
  if NOT(pos(ssGlobal,var)) then
    AVars(var):=val;
  end;

  ssCurSys(4,0):=var;
  ssCurSys(5,0):="";
  ssCurSys(6,0):=val;
  ssCurSys(7,0):=0;
end;


// check if equations are valid
// assumes equations have been saved
export ssCheckEqns()
begin
  local n:=size(ssCurSys(2));
  local eqn:="";
  local j,k,v;
  local varlist:={};
  local var:="";

  if n then

    // check each equation
    for j from 1 to n do
      eqn:=ssCurSys(2,j);

      // do we have = symbol?
      if instring(eqn,"=") then
        iferr eqn:=string(E0,1,12); then eqn:=""; end;
        E0:=ssCurSys(2,j);
        varlist:=LNAME(E0);
        v:=size(varlist);
        if v then
          for k from 1 to v do
            var:=string(varlist(k));
            if (pos(ssCurSys(4),var) == 0) then
              if ssAutoVar then
                ssAutoAddVar(var);
              else
                ssErrorPage(ip((j-1)/10)+1,ssNonVar,ssCurSys(2,j)+"\n> "+varlist(k));
                iferr E0:=eqn; then end;
                return(0);
              end;
            end;
          end;
        else
          // constants only? not equation
          ssErrorPage(j,ssNotEqn,ssCurSys(2,j));
          iferr E0:=eqn; then end;
          return(0);       
        end;
        iferr E0:=eqn; then end;
        
      else
        // no equal symbol; report page
        ssErrorPage(j,ssNotEqn,eqn);
        return(0);
      end;
    end;

    // all eqns seem ok
    return(1);

  else
    // no equations to check
    return(-1);
  end;
end;


// save current system
// no null names but does not check
// validity of names beyond whether
// AVars will accept such a name
ssSaveSystem()
begin
  ssSaveEqns();

  local sysname:=ssCurSys(1);
  local eqns:=ssCurSys(2);
  local vars:=ssCurSys(4);
  local eqn,var,n;
  local run:=1;

  // bug in input when choose lists are empty
  if (size(eqns) < 1) then eqns:={ ssVNoEqns }; end;
  if (size(vars) < 1) then vars:={ ssVNoVars }; end;


  while run do
    if input(
      {
        { sysname, [2], { 15, 75, 0 } },
        { eqn, eqns, { 15, 75, 1 } },
        { var, vars, { 15, 75, 2 } }
      },
      ssTSaveSys,
      { ssLName,    ssLEqns, ssLVars },
      { ssHSysName, ssHEqns, ssHVars }
    )
    then

      if size(sysname) then
        // maybe later add check for dupes
        n:=pos(ssSysTitles, sysname);

        if ssCheckEqns() then
          ssCurSys(1):=sysname;
          ssCurLib(ssCurSysIndex):=ssCurSys;
          AFiles(ssLibFN):=ssCurLib;
          if (ssCurSysIndex == 0) then ssCurSysIndex:=size(ssCurLib); end;
          ssLibSize:=size(ssCurLib);
        else
          ssWarn(ssNotSaved);
        end;
        run:=0;

      else
        msgbox(ssNullName);
      end;

    else
      // user canceled
      ssWarn(ssNotSaved);
      run:=0;
    end;

  end; // end while;
end;

// create a system
ssNewSystem()
begin
  ssCurSys:={ssTNewSys, {}, {}, {}, {}, {}, {} };
  ssCurSysIndex:=0;
  ssUpdateSys();
  ssInitVars();
  ssLoadSymb();
  startview(0,1);
end;


ssChooseSys()
begin
  local n,t;
  local run:=1;
  local systitles:={};

  n:=ssCurSysIndex;

  iferr
    ssCurLib:=AFiles(ssLibFN);
  then
    ssWarn(ssNoLibDat);
    ssCurLib:={};
    AFiles(ssLibFN):=ssCurLib;
  end;
  ssLibSize:=size(ssCurLib);

  if ssLibSize then
    systitles:=makelist(ssCurLib(n,1),n,1,ssLibSize);
  else
    systitles:={};
  end;
  ssSysTitles:=systitles;
  systitles(0):=ssTNewSys;

  t:=choose(n, ssTChooseSys, systitles);
  if t then
    if ( n <= ssLibSize ) then
      ssSetSystem(n);
    else
      ssNewSystem();
    end;
    return(1);
  end;
  return(0);
end;


// UI for deleting systems
export ssDeleteSys()
begin
  local cmd:="input({";
  local systitles;
  local j;
  local selected:={};
  local lib:={};

  if ssLibSize then
    selected:=makelist(0,j,1,ssLibSize);
    systitles:=makelist(ssCurLib(j,1),j,1,ssLibSize);
    for j from 1 to ssLibSize do
      cmd:=cmd + "{selected(" + string(j,1,0) + "),0,{94,5," + string(j-1,1,0) + "}}";
      if (j<ssLibSize) then cmd:=cmd + ",\n"; end;
    end;
    cmd:=cmd + "},\n" + string(ssTDeleteSys) + ",\n{";
    for j from 1 to ssLibSize do
      cmd:=cmd + string(ssCurLib(j,1));
      if (j<ssLibSize) then cmd:=cmd + ",\n"; end;
    end;
    cmd:=cmd + "},\n{";
    for j from 1 to ssLibSize do
      cmd:=cmd + string(ssHMarkDel);
      if (j<ssLibSize) then cmd:=cmd + ",\n"; end;
    end;
    cmd:=cmd + "},\n" + string(selected,1,0) + ",\n" + string(selected,1,0) + ")";
    Notes("eqlib_debug.log"):=cmd;
    j:=expr(cmd);

    if j then

      // current system deleted? if so set index to 0
      if ssCurSysIndex then
        if selected(ssCurSysIndex) then
          ssCurSysIndex:=0;
        end;
      end;

      // rebuild library and adjust index
      for j from 1 to ssLibSize do
        if NOT(selected(j)) then
          lib(0):=ssCurLib(j);
          if (j < ssCurSysIndex) then
            ssCurSysIndex:= ssCurSysIndex - 1;
          end;
        end;
      end;
      ssCurLib:=lib;
      AFiles(ssLibFN):=ssCurLib;

    end;

  else
    ssWarn(ssNoSystems);
  end;
end;


// variable browser
export ssVarBrowser()
begin
  local vals:=ssCurSys(6);
  local cons:=ssCurSys(7);
  local cmd:="input({";
  local n:=size(vals);
  local j,k;
  local tmp:="";
  local xloc:="15,75";

  if (n<1) then
    msgbox(ssNoVars);
    return(0);
  end;

  k:=1;
  for j from 1 to n do
    k:=max(k,size(ssCurSys(4,j)));
  end;

  if (k > 3) then xloc:="30,60"; end;

  // build input( ... )

  // { variable, type, { pos } }
  for j from 1 to n do
    k:=pos(ssGlobal, ssCurSys(4,j));
    if (k<28) then tmp:="0"; else tmp:="0,3"; end;
    
    cmd:=cmd + "{vals(" + string(j,1,0) + "),[" + tmp + "],{" + xloc + "," + string(j-1,1,0) + "}}, ";
    cmd:=cmd + "{cons(" + string(j,1,0) + "),0,{94,5," + string(j-1,1,0) + "}}";
    if ((n>1) AND (j<n)) then cmd:=cmd + ",\n"; end;
  end;

  // add title
  cmd:=cmd + "},\n" + string(ssCurSys(1) + " Variables") + ",\n{";

  // add labels
  for j from 1 to n do
    cmd:=cmd + string(ssCurSys(4,j) + "=") + ", " + string("");
    if ((n>1) AND (j<n)) then cmd:=cmd + ",\n"; end;
  end;
  cmd:=cmd + "},\n{";

  // add help
  for j from 1 to n do
    cmd:=cmd + string(ssCurSys(4,j) + ": " + ssCurSys(5,j)) + "," + string(ssHMakeConst);
    if ((n>1) AND (j<n)) then cmd:=cmd + ",\n"; end;
  end;
  cmd:=cmd + "},\n";

  // add init and reset vals
  tmp:="{";
  for j from 1 to n do
    tmp:=tmp + vals(j) + "," + cons(j);
    if ((n>1) AND (j<n)) then tmp:=tmp + ","; end;
  end;
  tmp:=tmp + "}";
  cmd:=cmd + tmp + ",\n" + tmp + ")";

  Notes("eqlib_debug.log"):=cmd;
  j:=expr(cmd);

  if j then
    ssCurSys(6):=vals;
    ssCurSys(7):=cons;

    for j from 1 to n do
      if pos(ssGlobal,ssCurSys(4,j)) then
        expr(ssCurSys(4,j) + ":=" + vals(j));
      else
        AVars(ssCurSys(4,j)):=vals(j);
      end;
    end;
  end;

  return(j);

end;



// build F, F(var1, var2, ...) and Jacobian
// assumes all equations and variables have been
// validated; also assumes equations are saved and
// that selectd equations form a non-null list
// saves F as ssF (CAS variable)
export casMakeF()
begin
  local cmd:="['";
  local cmd2:="[";
  local j,n;
  local varlist:={};
  local eqnlist:={};
  local k:=size(ssCurSys(4));
  local vstr1:="";
  local vstr2:="";

  n:=size(ssCurSys(2));
  for j from 1 to n do
    if ssCurSys(3,j) then
      eqnlist(0):=ssCurSys(2,j);
    end;
  end;

  for j from 1 to k do
    if NOT(ssCurSys(7,j)) then
      varlist(0):=j;
    end;
  end;
  ssSolveVars:=varlist;

  n:=k-1;
  for j from 1 to n do
    vstr1:=vstr1 + ssCurSys(4,j) + ",";
  end;
  vstr1:=vstr1 + ssCurSys(4,k);

  n:=size(varlist)-1;
  for j from 1 to n do
    vstr2:=vstr2 + ssCurSys(4,varlist(j)) + ",";
  end;
  vstr2:=vstr2 + ssCurSys(4,varlist(n+1));

  // ssF:=( vars ) -> [ f1, f2, ... ]
  n:=size(eqnlist)-1;
  for j from 1 to n do
    cmd:=cmd + eqnlist(j) + ")','";
    cmd2:=cmd2 + eqnlist(j) + "),";
  end;
  cmd:=cmd + eqnlist(n+1) + ")']";
  cmd:=replace(cmd,"=","-(");
  CAS(EVAL("ssF:=" + cmd));
  cmd2:=cmd2 + eqnlist(n+1) + ")]";
  cmd2:=replace(cmd2,"=","-(");
  CAS(EVAL("ssFx:=(" + vstr1 + ")->(" + cmd2 + ")"));

  cmd:="ssJ:=transpose(diff(" + cmd + ",[" + vstr2 + "]))";
  CAS(EVAL(cmd));
end;


// pseudo inverse using SVD
export ssPInv()
begin
  local l:=CAS("svd(ssJ)");
  local r:=size(l(2)); r:=r(1);
  local j,n,ns,s;
  local tol:=1E-30;
  local d:=[[0]];

  // built-in SVD fails on some singular matrices
  // failed SVD matrices returned as list of lists
  if ((type(l(1)) == 6) OR (type(l(3)) == 6)) then
    d:=CAS("pinv(ssJ,1e-30)");
    return(d);
  end;

  // sort sing. vals from largest to smallest
  // and adjust col/row accordingly
  for n from r downto 2 do
    ns:=1;
    s:=l(2,1);
    for j from 2 to n do
      if (l(2,j) <= s) then
        ns:=j; s:=l(2,j);
      end;
    end;
    if (ns < n) then
      l(2,ns):=l(2,n);
      l(2,n):=s;
      d:=l(1);
      l(1):=swapcol(d,ns,n);
      d:=l(3);
      l(3):=swapcol(d,ns,n);
    end;
  end;

  for n from 1 to r do
    s:=l(2,n);
    if (s >= tol) then
      l(2,n):=1/s;
    else
      l(2,n):=0;
    end;
  end;

  j:=size(l(3)); n:=j(1);
  j:=size(l(1)); r:=j(1);
  d:=makemat(ifte(I==J,l(2,I),0),n,r);

  // make sure to use CAS.trn for Hermetian
  return(l(3)*d*CAS.trn(EVAL(l(1))));
end;



// Newton's Method
// dx = - J^(-1) * F
// x[n+1] = x[n] + dx
export casNewton()
begin
  local j,k,nsv,nv;
  local var;
  local tolx:=1.0E-8;
  local tolF:=1.0E-10; // exit if F.F < tolF
  local nsF=2*tolF; // initial squared norm of F
  local vecNewt:={};
//  L1:={}; // for debugging only
  local pinv:=[[0]];
  local dx:=2*tolx;

  nv:=size(ssCurSys(4));
  nsv:=size(ssSolveVars);
  for j from 1 to nsv do
    vecNewt(0):=ssCurSys(6,ssSolveVars(j));
//    L1(0):=ssCurSys(6,ssSolveVars(j));
  end;
  
  for k from 1 to 100 do
    nsF:=CAS("approx(sqrt(ssF*ssF))");
    if (type(nsF) == 2) then
      msgbox(ssFUndef);
      break;
    end;
    vecNewt(0):=nsF;
    if (k == 1) then 
      vecNewt(0):=0;
    else
      vecNewt(0):=abs(dx);
    end;
    if (nsF < tolF) then break; end;
    if (abs(dx) < tolx) then break; end;

    iferr
      pinv:=ssPInv();
      dx:=-pinv*CAS.transpose("approx(ssF)");
    then
      msgbox(ssBadJac);
      break;
    end;

    for j from 1 to nsv do
      var:=ssCurSys(4,ssSolveVars(j));
      if pos(ssGlobal,var) then
        expr(var + ":=" + var + "+" + dx(j,1));
      else
        AVars(var):=AVars(var) + dx(j,1);
      end;
      ssCurSys(6,ssSolveVars(j)):=expr(var);
    end; // for j

    for j from 1 to nsv do
      vecNewt(0):=ssCurSys(6,ssSolveVars(j));
//      L1(0):=ssCurSys(6,ssSolveVars(j));
    end;
  end; // for k

  return(CAS.list2mat(EVAL(vecNewt),EVAL(nsv+2)));

end;


START()
begin
//  startview(-1,1);
  ssInitSSApp();
  ssInitVars();
  startview(0,1);
end;


// Equations Editor
Symb()
begin
  if ssInit then
    ssSaveEqns();
  else
    ssInitSSApp();
    ssInitVars();
  end;
  ssLoadSymb();
  startview(0,1);
end;


// System of Equations browswer
Plot()
begin
  if ssInit then
    ssSaveSystem();    
  else
    ssInitSSApp();
    ssInitVars();
  end;
  ssChooseSys();
  startview(0,1);
end;


// run solver
Num()
begin
  local matNewt:=[[0]];
  local j;
  local tlist:={};
  local labels:={};
  if ssInit then
    ssSaveEqns();
    if NOT(ssCheckEqns()) then
      return(0);
    end;
    if NOT(ssHaveEqns()) then
      return(0); // exit
    end;

    if ssVarBrowser() then
      casMakeF();
      matNewt:=casNewton();
      tlist(0):=ssCurSys(1);
      tlist(0):={};
      for j from 1 to size(ssSolveVars) do
        labels(0):=ssCurSys(4,ssSolveVars(j));
      end;
      labels(0):="ǁFǁ";
      labels(0):="ǁΔxǁ";
      tlist(0):=labels;
      editmat(matNewt,tlist);
    end;
  else
    ssInitSSApp();
    ssInitVars();
  end;
end;


// Set up initial guess 
NumSetup()
begin
  if ssInit then ssSaveEqns(); end;
  ssVarBrowser();
end;


// select equations page (View Menu)
view "Select Page", ssMenuEqnPage()
begin
  ssSelectEqPage();
  startview(0,1);
end;


// UI for creating new variables (View Menu)
view "Add/Edit Variable", ssMenuNewVar()
begin
  ssCreateVar();
end; 

// reload saved equations (View Menu)
view "Restore Equations", ssMenuLoadEqns()
begin
  if ssCurSysIndex then
    ssCurSys:=ssCurLib(ssCurSysIndex);
    ssUpdateSys();
    ssInitVars();
  end;
  ssLoadSymb();
  startview(0,1);
end;


// save current system (View Menu)
view "Save System", ssMenuSaveSys()
begin
  ssSaveSystem();
end;

// create new system (View Menu)
view "New System", ssMenuNewSys()
begin
  ssSaveSystem();
  ssNewSystem();
end;

// delete existing system (View Menu)
view "Delete System", ssMenuDelSys()
begin
  ssSaveSystem();
  ssDeleteSys();
end;

// select system (View Menu)
view "Select System", ssMenuSetSys()
begin
  ssSaveSystem();
  if ssChooseSys() then startview(0,1); end;
end;

Original first post below -- info is out of date due to revisions
-------------------------------------------------------------------------------------------
Here's a brief outline of my ideas.

1) Each category of equations will be a list of lists of the following form:

Code:

{
  "Title",
  { "equation1", "equation2", ..., "equationn" },
  { "var1", "var2", ..., "varn" },
  { "vardescription1", "vardescription2", ..., "vardescriptionn" }
}

Perhaps later on we can include graphics by means of lists of point definitions, line definitions, and polygon definitions. Ideally, the user would not have to specify what the variable list is, and the app will be smart enough to parse through each equation string to identify what the variables are.

2) A choose screen will allow users to select which category of equations to work with

3) Once a category is selected, all the equations will be shown in an INPUT() window with check boxes next to them. The user then checks the equations to use.

4) The app then creates another INPUT() window for variables and their initial values. Each variable may also be checked to specify whether the quantity should be treated as a variable that may change during the solve process, or if it is a constant (much like in the SOLVESYS library for the HP48 and HP49).

5) The equations are then modified to insert constants where appropriate, and then call fsolve() to solve the system.

This is the general scheme. What would also be nice is an interface for users to add in new categories of equations. Anyway, this is something for a later update. For now, I'll start working on building the initial interface with the assumption that categories of equations already exist. This is where, if you want to contribute, you can help by posting the categories as shown in the format above.

I'll start with a simple example:

Code:

{
  "Linear Motion",
  { 
    "x=x0+v0*t+1/2*a*t^2",
    "x=x0+v*t-1/2*a*t^2",
    "x=x0+1/2*(v0+v)*t",
    "v=v0+a*t"
  },
  { "x", "x0", "v0", "v", "t", "a" },
  {
    "final position",
    "initial position",
    "initial velocity",
    "final velocity",
    "time",
    "acceleration"
  }
}

[ArielPalazzesi]

Hi Han!

I like the way you've addressed the development of the library. I'll see your progress closely .

A hug .

[Han]

Here's a very basic skeleton. This version does not handle systems that have free variables. Ideally, in the case of free variables, then the initial guesses should be used to establish a solution. For example, if (x,y,z) = (x,2x-1,3x+1), then the initial guess of (4,5,6) should produce the solution (4,7,13) using x=4. This will be implemented over time. For now, what you get is basically a proof of concept.

Code:

export eqldata;
export varlist;
export eqlist;
export eql_results;
export inputstr;

eql_init();
eql_cmenu();
eql_setup();
eql_dovars();
eql_dosolve();

export eqlib()
begin
  eql_init();
  local c:=eql_cmenu();
  local tmp;

  if c then
    tmp:=eql_setup(c);
  end;

  if tmp then
    tmp:=eql_dovars(c,1);
  end;

  if tmp then
    repeat
      eql_dosolve(c);
      tmp:=eql_dovars(c,0);
    until tmp==0;
  end;

end;

eql_init()
begin
  local tmp;

  iferr tmp:=Notes("EqLib") then
    msgbox("No equation library data!");
  else
    eqldata:=expr(tmp);
  end;
end;

eql_cmenu()
begin
  local c;
  local n:=size(eqldata);
  local categories:=makelist(eqldata(X,1),X,1,n);

  choose(c,"Equation Library", categories);

  return(c);
end;

eql_setup(c)
begin
  local eqns:=eqldata(c);
  local tmp:="input({";
  local n:=size(eqns(2));
  local i;
  eqlist:=makelist(1,X,1,n);

  for i from 1 to n do
    tmp:=tmp+"{eqlist(" + string(i,1,12) + "),0,{94,5," + string(i-1,1,12) + "}}";
    if (n>1) AND (i<n) then tmp:=tmp + ","; end;
  end;
  tmp:=tmp + "}," + string(eqns(1)) + ",{";

  for i from 1 to n do
    tmp:=tmp + string(eqns(2,i));
    if (n>1) AND (i<n) then tmp:=tmp + ","; end;
  end;
  tmp:=tmp + "},{";

  for i from 1 to n do
    tmp:=tmp + string("Use this equation?");
    if (n>1) AND (i<n) then tmp:=tmp + ","; end;
  end;
  tmp:=tmp + "})";
  n:=expr(tmp);
  return(n);
end;

eql_dovars(c,f)
begin
  local eqns:=eqldata(c);
  local tmp:="";
  local n:=size(eqns(3));
  local i; local t;
  if f then
    varlist:=makelist(0,X,1,n);
  end;

  for i from 1 to n do
    tmp:="type(" + eqns(3,i) + ")";
    if CAS(tmp)==6 then tmp:=eqns(3,i) + ":=0"; CAS(tmp); end;
  end;

  tmp:="input({";

  for i from 1 to n do
    tmp:=tmp + "{" + eqns(3,i) + ",[0],{30,30," + string(i-1,1,12) + "}},";
    tmp:=tmp + "{varlist(" + string(i,1,12) + "),0,{94,5," + string(i-1,1,12) + "}}";
    if (n>1) AND (i<n) then tmp:=tmp + ","; end;
  end;
  tmp:=tmp + "}," + string(eqns(1) + " Variables") + ",{";

  for i from 1 to n do
    tmp:=tmp + string(eqns(3,i) + "=") + "," + string("");
    if (n>1) AND (i<n) then tmp:=tmp + ","; end;
  end;
  tmp:=tmp + "},{";

  for i from 1 to n do
    tmp:=tmp + string("Enter the value for the " + eqns(4,i)) + ",";
    tmp:=tmp + string("Make this variable a constant?");
    if (n>1) AND (i<n) then tmp:=tmp + ","; end;
  end;
  tmp:=tmp + "},";

  local tmp2:="{";

  for i from 1 to n do
    tmp2:=tmp2 + CAS(eqns(3,i)) + "," + varlist(i);
    if (n>1) AND (i<n) then tmp2:=tmp2 + ","; end;
  end;
  tmp2:=tmp2 + "}";

  tmp:=tmp + tmp2 + "," + tmp2 + ")";
  t:=expr(tmp);

  if t then
    for i from 1 to n do
      if varlist(i)==0 then
        tmp:="purge(" + eqns(3,i) + ")";
        CAS(tmp);
      end;
    end;
  end;
  
  return(t);
end;

eql_dosolve(c)
begin
  local eqns:=eqldata(c);
  local n:=size(eqns(2));
  local i;
  local tmp:="fsolve([";
  local k:=ΣLIST(eqlist);

  for i from 1 to n do
    if eqlist(i) then
      k:=k-1;
      tmp:=tmp + eqns(2,i);
      if k then tmp:=tmp + ","; end;
    end;
  end;
  tmp:=tmp + "],[";

  n:=size(eqns(3));
  k:=n-ΣLIST(varlist);

  for i from 1 to n do
    if varlist(i)==0 then
      k:=k-1;
      tmp:=tmp + eqns(3,i);
      if k then tmp:=tmp + ","; end;
    end;
  end;
  tmp:="eql_results:=" + tmp + "])";
  inputstr:=tmp;
  CAS(tmp);
  local l:=size(eql_results);

  if type(eql_results)==6 then msgbox("No solution or infinite solutions; not implemented"); end;

// if multiple solutions; just get the first one; fix this later to choose "best" solution
  if type(eql_results)==4 then
    k:=1;
    for i from 1 to n do
      if varlist(i)==0 then
        tmp:=eqns(3,i) + ":=" + eql_results(1,k);
        CAS(tmp);
        k:=k+1;
      end;
    end;
  end;

end;

You must first create a Note named "EqLib" whose contents is:

Code:

{
  {
    "Linear Motion",
    { 
      "x=x0+v0*t+1/2*a*t^2",
      "x=x0+v*t-1/2*a*t^2",
      "x=x0+1/2*(v0+v)*t",
      "v=v0+a*t"
    },
    { "x", "x0", "v0", "v", "t", "a" },
    {
      "final position",
      "initial position",
      "initial velocity",
      "final velocity",
      "time",
      "acceleration"
    }
  }
}

Run eqlib() and choose "Linear Motion" (the only option). The next screen allows you to pick which equations to use. Leave all selected. The screen after that allows you to set up the known values. Those which should be treated as constants can be checked so that they are not considered a variable in the system. Press OK. If a solution exists, the values are updated. Press cancel to quit. I chose x=1, x0=5, t=3, and a=9.8. This results in v0=-16.033333... and v=13.366666...

Again, this is just a proof of concept. There are still a number of things I have to figure out about fsolve. Even when a system is solvable, if we give fsolve initial conditions, we sometimes get "Invalid dimension" errors. It appears that if fsolve is given more equations than the number of variables, it claims invalid dimension even though that may not be the case (e.g. some of the equations are equivalent).

An alternative to fsolve would be to implement our own system solver using something similar to that in SOLVESYS for the HP48/HP49. A brief glance at the source code suggests a Newton-like algorithm using the Jacobian.

[Giancarlo]

Hello Han,
thank you so much for working on an equation library. I tried to use it in the EMU. Everything went fine but i receive the "not supported" msgbox for any value assigned to the variables.

Apart from this I would like to understand why you chose to store the equations in a note and not in a simple list. Which advantages you see in doing like this?

Second question, why do you prefer to use fsolve() instead of programming the solve app? just to manage the overall process?

Thanks

Giancarlo

[Han]

(01-10-2015 04:12 PM)Giancarlo Wrote:  Hello Han,
thank you so much for working on an equation library. I tried to use it in the EMU. Everything went fine but i receive the "not supported" msgbox for any value assigned to the variables.

Apart from this I would like to understand why you chose to store the equations in a note and not in a simple list. Which advantages you see in doing like this?

Second question, why do you prefer to use fsolve() instead of programming the solve app? just to manage the overall process?

Thanks

Giancarlo

If the system and initial guesses result in anything other than a single solution, then the generic "not supported" message will appear. The current skeleton of an engine does not handle any case beyond the single-solution case (for the moment; that is not the plan for the future, of course). If you simply entered initial guesses without specifying any of the "variables" to be constants, then there would be free variables (since we have more variables than equations), resulting in a system with infinite solutions.

As for why I stored the data in notes -- it's to prevent the erasure of data when one opens the source code of the program. Any variable associated with the program will be reset. Now, we could embed the data into the source code, but as you know, pressing [Shift][Esc] erases the entire program without any prompts. I don't know how many times I have cursed the virtual calculator or the real calculator because I didn't notice that the [Shift] operator was already active and accidentally erased my entire program. Saving it in the notes means we avoid this issue. Since the plan is to also implement an interface to add equations and manage the library, there is no real reason for any user to manually edit the corresponding notes (since it too has the same issue as programs with the [Shift][Esc]).

[akmon]

Hello Han. Your Project is very interesting, and believe me, I was thinking of programming something similar to Solvesys. In my case using a process to extract the variables of the equations, avoiding to make a variable list, excepting if you deliberately insert a list explaining its purpose, like your example. So it´s a good idea to extract them automatically, as you wrote. Don´t know if there is a command like TVARS on HP 49 to get the variables of a formula.

By the way, I cannot find the solutions on you examples in my emulator. After inserting values a Box appears saying: no solutions or infinite solutions, and then time value returns to 0. I´ll try on my physical calculator.

Another suggestion, the possibility to go back to the formula menu from the variable menu.

And finally, fsolve Works better than solve App, let´s see if it behaves well with more examples. With this update (6975), I managed to solve a system of equations using fsolve where solve App couldn´t.

[Han]

To obtain the answers I got, you must specify that x, x0, t, and a are constants (check box) next to the variable.

Edit: I see the issue; there is a priority conflict between local and CAS variables. This should be fixed relatively soon.

[akmon]

I did it, but no luck.

[Han]

Try this slight update:

Code:

export eql_data;
export eql_varlist;
export eql_eqnslist;
export eql_initvals;
export eql_results;

// this is for debugging
export inputstr;

eql_init();
eql_cmenu();
eql_setup();
eql_dovars();
eql_dosolve();

export eqlib()
begin
  eql_init();
  local c:=eql_cmenu();
  local tmp;

  if c then
    tmp:=eql_setup(c);
  end;

  if tmp then
    tmp:=eql_dovars(c,1);
  end;

  if tmp then
    repeat
      eql_dosolve(c);
      tmp:=eql_dovars(c,0);
    until tmp==0;
  end;

end;

eql_init()
begin
  local tmp;

  iferr tmp:=Notes("EqLib") then
    msgbox("No equation library data!");
  else
    eql_data:=expr(tmp);
  end;
end;

eql_cmenu()
begin
  local c;
  local n:=size(eql_data);
  local categories:=makelist(eql_data(X,1),X,1,n);

  choose(c,"Equation Library", categories);

  return(c);
end;

eql_setup(c)
begin
  local eqns:=eql_data(c);
  local tmp:="input({";
  local n:=size(eqns(2));
  local i;
  eql_eqnslist:=makelist(1,X,1,n);

  for i from 1 to n do
    tmp:=tmp+"{eql_eqnslist(" + string(i,1,12) + "),0,{94,5," + string(i-1,1,12) + "}}";
    if (n>1) AND (i<n) then tmp:=tmp + ","; end;
  end;
  tmp:=tmp + "}," + string(eqns(1)) + ",{";

  for i from 1 to n do
    tmp:=tmp + string(eqns(2,i));
    if (n>1) AND (i<n) then tmp:=tmp + ","; end;
  end;
  tmp:=tmp + "},{";

  for i from 1 to n do
    tmp:=tmp + string("Use this equation?");
    if (n>1) AND (i<n) then tmp:=tmp + ","; end;
  end;
  tmp:=tmp + "})";
  n:=expr(tmp);
  return(n);
end;

eql_dovars(c,f)
begin
  local eqns:=eql_data(c);
  local tmp:="";
  local n:=size(eqns(3));
  local i; local t;
  if f then
    eql_varlist:=makelist(0,X,1,n);
    eql_initvals:=makelist(0,X,1,n);
  end;

  tmp:="input({";

  for i from 1 to n do
    tmp:=tmp + "{eql_initvals(" + string(i,1,12) + "),[0],{30,30," + string(i-1,1,12) + "}},";
    tmp:=tmp + "{eql_varlist(" + string(i,1,12) + "),0,{94,5," + string(i-1,1,12) + "}}";
    if (n>1) AND (i<n) then tmp:=tmp + ","; end;
  end;
  tmp:=tmp + "}," + string(eqns(1) + " Variables") + ",{";

  for i from 1 to n do
    tmp:=tmp + string(eqns(3,i) + "=") + "," + string("");
    if (n>1) AND (i<n) then tmp:=tmp + ","; end;
  end;
  tmp:=tmp + "},{";

  for i from 1 to n do
    tmp:=tmp + string("Enter the value for the " + eqns(4,i)) + ",";
    tmp:=tmp + string("Make this variable a constant?");
    if (n>1) AND (i<n) then tmp:=tmp + ","; end;
  end;
  tmp:=tmp + "},";

  local tmp2:="{";

  for i from 1 to n do
    tmp2:=tmp2 + eql_initvals(i) + "," + eql_varlist(i);
    if (n>1) AND (i<n) then tmp2:=tmp2 + ","; end;
  end;
  tmp2:=tmp2 + "}";

  tmp:=tmp + tmp2 + "," + tmp2 + ")";
  t:=expr(tmp);

  if t then
    for i from 1 to n do
      if eql_varlist(i) then
        tmp:=eqns(3,i) + ":=" + eql_initvals(i);
      else
        tmp:="purge(" + eqns(3,i) + ")";
      end;
      CAS(tmp);
    end;
  end;
  
  return(t);
end;

eql_dosolve(c)
begin
  local eqns:=eql_data(c);
  local n:=size(eqns(2));
  local i;
  local tmp:="fsolve([";
  local k:=ΣLIST(eql_eqnslist);

  for i from 1 to n do
    if eql_eqnslist(i) then
      k:=k-1;
      tmp:=tmp + eqns(2,i);
      if k then tmp:=tmp + ","; end;
    end;
  end;
  tmp:=tmp + "],[";

  n:=size(eqns(3));
  k:=n-ΣLIST(eql_varlist);

  for i from 1 to n do
    if eql_varlist(i)==0 then
      k:=k-1;
      tmp:=tmp + eqns(3,i);
      if k then tmp:=tmp + ","; end;
    end;
  end;
  tmp:="eql_results:=" + tmp + "])";
  inputstr:=tmp;
  CAS(tmp);
  local l:=size(eql_results);

  if type(eql_results)==6 then msgbox("No solution or infinite solutions; not implemented"); end;

// if multiple solutions; just get the first one; fix this later to choose "best" solution
  if type(eql_results)==4 then
    k:=1;
    for i from 1 to n do
      if eql_varlist(i)==0 then
        eql_initvals(i):=eql_results(1,k);
        k:=k+1;
      end;
    end;
  end;

end;

The issue was that the source code used the local variable 't' which also happens to be a variable that we declare as a CAS variable via indirection. Fortunately (or unfortunately depending how you view the separation of Home and CAS), the CAS variables seem to be completely separate. Therefore, we can represent the value of the variables as entries within a list we call eql_initvals. Then, by declaring variables (or purging) within the CAS and setting them equal to eql_initvals(i), we do not run into any collision between local and non-local variables.

This scheme can be broken if for some reason we decide to create equations whose variables are also of the same name as those used within our source code! It is my hope, however, that the variable names I used are unique enough that they never appear in any useful equation that people use in real work :-)

[akmon]

Now it Works fine. Thank you for the new source code.

Another suggestion.
Is it posible to add an equation directly on the program?

I think I´ve left some clues to understand I am a big fan of SOLVESYS library, my mind is forcing to suggest your program to be as similar as posible to that one. My apologies.
The purpose is to have a program that can read a note file with all equations (eqlib), or simply make a fast calculus inserting directly one or several equations for an specific problem (solvesys). All in one code. fsolve is better than solve App for numeric solution (a high percentaje of complex equations problems you just need the final numeric value).

[Han]

(01-10-2015 10:49 PM)akmon Wrote:  Now it Works fine. Thank you for the new source code.

Another suggestion.
Is it posible to add an equation directly on the program?

I think I´ve left some clues to understand I am a big fan of SOLVESYS library, my mind is forcing to suggest your program to be as similar as posible to that one. My apologies.
The purpose is to have a program that can read a note file with all equations (eqlib), or simply make a fast calculus inserting directly one or several equations for an specific problem (solvesys). All in one code. fsolve is better than solve App for numeric solution (a high percentaje of complex equations problems you just need the final numeric value).

Yeah, the idea was to basically make a hybrid of SOLVESYS and the EQLIB from the HP48 days.

[acmeanvil]

Hi Han,
That is an excellent prototype. I will follow progress closely.

And you are correct, it would be great to be able to incorporate images in with the choose box...

[Han]

The solver portion is done and behaves much like the SOLVESYS solver. The slight difference is that the solver falls back to fsolve in the case where there are more variables than there are equations (usually infinite solutions). In this case, the non-constant initial values are ignored. If there are more equations than variables, then a least squares approach will be used.

What I am working on now is the interface. Sadly, the CHOOSE() and MSGBOX() commands are just not sufficient for what I need. Moreover, they keep switching the view back to the command line view rather than just drawing directly onto the current view.

[akmon]

Fantastic, I´m waiting for the code to test some linear and no linear system of equations.

[Helge Gabert]

Nice project, Han!

Do you plan on having check boxes for variables that are allowed to be complex as well? (The user might want the solution(s) to be restricted to reals in some cases, but at other times not).

[Han]

After some testing, I have run into a snag. Apparently the recent firmware changed how the command SVD outputs its results -- and the results themselves!

[Helge Gabert]

I noticed that for arbitrary n x m matrices, the supposedly rectangular vector (the second element returned in the SVD output list), has dimension n if n<m, or m if n>m, so you can't just call diag(outputlist(2)) in order to reconstruct.

Is this the problem?

The results, in general, seem to be OK, though.

[Han]

No, that isn't the issue. Try this matrix:

SVD([[3,4.5,4.5],[4.5,12.25,1.25],[4.5,1.25,12.25]]);

Just copy that line and paste it into the command line of the emulator. I get some "noise" (complex-valued matrices and vectors); see attached screenshot. The bug is also present in hardware (latest firmware)

[Helge Gabert]

Yes I see what you are saying - the matrix is symmetric, and Wolfram alpha gives a purely real solution. Even with complex checked off in Home settings, and use i checked off in CAS settings it returns the complex elements in the U matrix.

But . . . presume the output is stored in L1, then, in CAS, L1(1)*diag(L1(2))*TRN(L1(3)) returns the original matrix - - so couldn't one argue that the complex SVD output is also correct?

[Han]

(01-14-2015 09:03 PM)Helge Gabert Wrote:  Yes I see what you are saying - the matrix is symmetric, and Wolfram alpha gives a purely real solution. Even with complex checked off in Home settings, and use i checked off in CAS settings it returns the complex elements in the U matrix.

But . . . presume the output is stored in L1, then, in CAS, L1(1)*diag(L1(2))*TRN(L1(3)) returns the original matrix - - so couldn't one argue that the complex SVD output is also correct?

Interesting... I didn't even bother to check that particular calculation. What I don't understand is why there is even a complex component when the characteristic polynomial has real roots.