c++: Solving dynamic system of nonlibear equations

Question

I'm trying to code model for calculating input rpms for two servo drives. The model consists of four equations with eight variables.

mean speed: v = (vd+vb)/ 2
slide to roll ratio: srr = (vd-vb)/ v
disc speed: vd = 2* pi* Rd* nd
ball speed: vb = 2* pi* Rb* nb

Inside application (QT gui), the model work like a table in which each cell represent variable. Some of the cells user fill with values and the rest is calculated based on the equations. The solution is done after each input - if user inserts "vd" and "vb", application solve "v" and "srr", then wait for the next input.

Now the model is solved by a bunch of ifs conditions defined for most likely possible combinations of inputs. This is however not a very elegant solution since the code is messy and hard to edit/extend.

Main problem/question is how to put these equations into a matrix and solve for every combination of inputs(input should be the definition of 4 variables).

//Class of the model
class Prediction : public QObject
{
    Q_OBJECT
public:
    Prediction(QObject* parent = nullptr);

    PhysicalQuantity constant(double value);

    //called to solve the model
    void evaluate();

signals:
    void callDataChanged();

private:
    const double pi = 3.141592653589793;
    int m_row;

    PhysicalQuantity m_meanSpeed  = PhysicalQuantity("v");
    PhysicalQuantity m_ballRPM = PhysicalQuantity("bRPM");
    PhysicalQuantity m_discRPM = PhysicalQuantity("dRPM");
    PhysicalQuantity m_ballSpeed = PhysicalQuantity("vb");
    PhysicalQuantity m_discSpeed = PhysicalQuantity("vd");
    PhysicalQuantity m_slideToRoll = PhysicalQuantity("srr");
    PhysicalQuantity m_discRadius = PhysicalQuantity("rd");
    PhysicalQuantity m_ballRadius = PhysicalQuantity("rb");
};

//current solution
void Prediction::evaluate()
{
    if(this->m_meanSpeed.isEntered() && this->m_slideToRoll.isEntered() && m_discRadius.isEntered() && m_ballRadius.isEntered())
    {
        m_discSpeed.computedValue(((m_slideToRoll * m_meanSpeed + 2 * m_meanSpeed)/m_slideToRoll).value());
        m_ballSpeed.computedValue((2*m_meanSpeed - m_discSpeed).value());
        m_discRPM.computedValue((m_discSpeed / (2 * pi * m_discRadius)).value());
        m_ballRPM.computedValue((m_ballSpeed / (2 * pi * m_ballRadius)).value());
    }
...
};

//class of the variables in the model
class PhysicalQuantity
{

public:
    PhysicalQuantity();
    PhysicalQuantity(QString id);
    PhysicalQuantity(double value);
    PhysicalQuantity(const PhysicalQuantity &obj);

    //read functions
    bool isEntered() const;
    bool isConstant() const;
    double value() const;
    double matInt() const;
    QString name() const;

    //write functions
    bool setEntered(bool isEnetered);
    QString name(QString id);
    double value(double value);
    double computedValue(double value);


    //Operators
    PhysicalQuantity operator+(const PhysicalQuantity &obj);
    PhysicalQuantity operator-(const PhysicalQuantity &obj);
    PhysicalQuantity operator/(const PhysicalQuantity &obj);
    PhysicalQuantity operator=(const PhysicalQuantity &obj);
    PhysicalQuantity operator-();

protected:
    bool m_isEntered;
    bool m_isConstant;
    bool m_isComputed;
    double m_value = 1234;
    double m_defaultVal;
    QString m_id;

Answer 1

Having

• v = (vd+vb)/ 2
• srr = (vd-vb)/ v
• vd = 2* pi* Rd* nd
• vb = 2* pi* Rb* nb

the possible formulas are ( I removed some requiring for nothing additional variable(s) like v = (vd-vb)/srr when there is already v = (vd+vb)/2 without needing srr) :

• vd = 2* pi* Rd* nd or vd = 2*v - vb

• vb = 2* pi* Rb* nb or vb = 2*v - vd

• v = (vd+vb)/2

• srr = (vd-vb)/v

• Rd = vd/(2* pi* nd)

• Rb = vb/(2* pi* nb)

• nd = vd / (2* pi* Rd)

• nb = vb/(2* pi* Rb)

so there are only 10 formulas, if you like programming by the data :

#include <string>
#include <vector>
#include <iostream>

const double pi = 3.141592653589793;

// index in the vector of variables
enum VARS { VD, VB, V, SSR, RD, RB, ND, NB, NVARS };

// manages a var
struct Var {
  std::string name;
  bool set;
  double value;
};

struct Formula {
  bool (*pf)(std::vector<Var> & vars); // the function to do the work
  VARS r; // the var set by the formula
  std::vector<VARS> used; // the vars used by the formula
};

// the functions doing computations, one per formula
// they are called only when the used vars are set
// return false if divide by 0

bool fvd1(std::vector<Var> & vars)
{
  vars[VD].value = 2*pi*vars[RD].value*vars[ND].value;
  return true;
}

bool fvd2(std::vector<Var> & vars)
{
  vars[VD].value = 2*vars[V].value - vars[VB].value;
  return true;
}

bool fvb1(std::vector<Var> & vars)
{
  vars[VB].value = 2*pi*vars[RB].value*vars[NB].value;
  return true;
}

bool fvb2(std::vector<Var> & vars)
{
  vars[VB].value = 2*vars[V].value - vars[VD].value;
  return true;
}

bool fv(std::vector<Var> & vars)
{
  vars[V].value = (vars[VD].value + vars[VB].value)/2;
  return true;
}

bool fssr(std::vector<Var> & vars)
{
  if (vars[V].value == 0)
    return false;

  vars[SSR].value = (vars[VD].value - vars[VB].value) / vars[V].value;
  return true;
}

bool fRd(std::vector<Var> & vars)
{
  if (vars[ND].value == 0)
    return false;

  vars[RD].value = vars[VD].value/(2 *pi * vars[ND].value);
  return true;
}

bool fRb(std::vector<Var> & vars)
{
  if (vars[NB].value == 0)
    return false;

  vars[RB].value = vars[VB].value/(2 *pi * vars[NB].value);
  return true;
}

bool fnd(std::vector<Var> & vars)
{
  if (vars[RD].value == 0)
    return false;

  vars[ND].value = vars[VD].value/(2 *pi * vars[RD].value);
  return true;
}

bool fnb(std::vector<Var> & vars)
{
  if (vars[RB].value == 0)
    return false;

  vars[NB].value = vars[VD].value/(2 *pi * vars[RB].value);
  return true;
}

const std::vector<Formula> Formulas = {
  { &fvd1, VD, { RD, ND } }, // to compute VD by fvd1 the needed vars are RD and ND
  { &fvd2, VD, { V, VB } },
  { &fvb1, VD, { RD, NB } },
  { &fvb2, VD, { V, VD } },
  { &fv, V, { VD, VB } },
  { &fssr, SSR, { VD, VB, V } },
  { &fRd, RD, { VD, ND } },
  { &fRb, RB, { VB, NB } },
  { &fnd, ND, { VD, RD } },
  { &fnb, NB, { VB, RB } }
};

// try to apply the formulas
// do in a loop in the case newly computed var(s) allows
// to set other(s) and that independently of the order 
// of the formulas
bool formulas(std::vector<Var> & vars, int & knownVars)
{
  bool modified;

  do {
    modified = false;

    for (const Formula & f : Formulas) {
      if (!vars[f.r].set) {
        bool ok = true;

        for (VARS v : f.used) {
          if (!vars[v].set) {
            ok = false;
            break;
          }
        }
        if (ok) {
          if (! (*f.pf)(vars))
            return false;
          vars[f.r].set = true;
          modified = true;
          knownVars += 1;
        }
      }
    }
  } while (modified);

  return true;
}

// to do in a terminal without graphic

int main()
{
  std::vector<Var> vars(NVARS);

  vars[VD] = { "vd", false, 0 };
  vars[VB] = { "vb", false, 0 };
  vars[V] = { "v", false, 0 };
  vars[SSR] = { "ssr", false, 0 };
  vars[RD] = { "Rd", false, 0 };
  vars[RB] = { "Rb", false, 0 };
  vars[ND] = { "nd", false, 0 };
  vars[NB] = { "nb", false, 0 };

  int knownVars = 0;

  do {
    std::cout << "enter var ( ";
    for (const Var & v : vars) {
      if (!v.set)
        std::cout << v.name << ' ';
    }
    std::cout << " ) ? ";

    std::string s;

    if (! (std::cin >> s)) {
      // EOF
      return -1;
    }

    std::vector<Var>::iterator it;

    for (it = vars.begin(); it != vars.end(); ++it) {
      if (!it->set && (it->name == s))
        break;
    }

    if (it == vars.end())
      std::cerr << "invalid name" << std::endl;
    else if (std::cout << "enter value ? ", ! (std::cin >> it->value)) {
      std::cerr << "invalid value" << std::endl;
      std::cin.clear();
      if (! (std::cin >> s)) {
        // EOF
        return -1;
      }
    }
    else {
      it->set = true;
      knownVars += 1;

      if (!formulas(vars, knownVars)) {
        std::cerr << "divide by 0, abort" << std::endl;
        return -1;
      }
    }
  } while (knownVars != NVARS);

  std::cout << "Done :";
  for (const Var & v : vars)
    std::cout << ' ' << v.name << '=' << v.value;
  std::cout << std::endl;

  return 0;
}

I do not use Qt but standard C++ string/vector/IO but that do not change the algorithm

Compilation and execution :

pi@raspberrypi:/tmp $ g++ -g -pedantic -Wextra -Wall f.cc
pi@raspberrypi:/tmp $ ./a.out
enter var ( vd vb v ssr Rd Rb nd nb  ) ? vd
enter value ? 1
enter var ( vb v ssr Rd Rb nd nb  ) ? nb
enter value ? 12
enter var ( vb v ssr Rd Rb nd  ) ? vb
enter value ? 2
enter var ( Rd nd  ) ? Rd
enter value ? 3
Done : vd=1 vb=2 v=1.5 ssr=-0.666667 Rd=3 Rb=0.0265258 nd=0.0530516 nb=12
pi@raspberrypi:/tmp $