How to change the equation into a code in an nth pendulum?

[Rafsing]

I want to simulate the nth level pendulum (double, triple and so on). I found the equations for the acceleration of the nth pendulum on this page: https://www.researchgate.net/publication/336868500_Equations_of_Motion_Formulation_of_a_Pendulum_Containing_N-point_Masses but I am having trouble changing it to code. My main problem is that I don’t know where some of the indexes came from. Here’s the equation I’d like to change into code:

This is my code:

 int sigma(int j, int k)
  {
    return int(j<=k);
  }

  int  phi(int j, int k)
  {
    return int(!(j==k));
  }


  void calculateAcceleration()
  {
    float dividend = 1, divider = 0;
    float[] copyAcceleration = new float[bob.size()];

    for (int i = 0; i <bob.size() - 1; i++)
    {
      copyAcceleration[i] = angleAcceleration.get(i);
    }

    for (int i = 0; i <bob.size() - 1; i++)
    {
      for (int j = 0; j <bob.size() - 1; j++)
      {
        dividend+= angles.get(j);
      }

      for (int k = 0; k <bob.size() - 1; k++)
      {
        divider += bob.get(k).mass * lenghts.get(bob.size() - 1) * lenghts.get(bob.size() - 1) * sigma(bob.size(), k);
      }


      copyAcceleration[i] = -dividend/divider;
    }

    for (int i = 0; i <bob.size() - 1; i++)
    {
      angleAcceleration.set(i, copyAcceleration[i]);
    }
  }

Here’s my whole project: https://pastebin.com/mAi1aNSL. Can someone help me with this?

[MiguelSanches]

Hi @Rafsing,

If I am not mistaken and didn’t commit any mistakes the following code should reproduce the equation

int n = 3;
float []m = {20, 20, 20};
float []l = {40, 40, 40};
float []angle = {PI/2, PI/2, PI/2};
float []vel = {0,0,0};
float []accel = {0,0,0};

void setup() {
  float g = 1;
  float den = 0;
  float num1 = 0;
  float num2 = 0;
  float num3 = 0;

  for (int j = 0; j < n; j++) {
    //denominator 
    for (int k = 0; k < n; k++) {
      den += m[k] * l[k] * l[k] * (j <= k? 1 : 0);
    }

    
    for (int k = 0; k < n; k++) {
      //first numerator
      num1 = g * l[j] * sin(angle[j]) * m[k] * (j <= k? 1 : 0);

      //second numerator
      float inner_sum = 0;
      // inner sum
      for (int q = k+1; q < n; q++) {
        inner_sum += m[q] * (j <= q? 1 : 0);
      }
      num2 = inner_sum * l[j] * l[k] * sin(angle[j] - angle[k]) * vel[j] * vel[k];
      
      //Third numerator
      //The inner sum is the same as in the num2
      num3 = inner_sum * l[j] * l[k] * (sin(angle[k] - angle[j]) * (vel[j] * vel[k]) * vel[k] + (j != k ? 1 : 0) * cos(angle[j] - angle[k])*accel[k]);
    }
    
    accel[j] = - (num1 + num2 + num3) / den;
  }
  printArray(accel);
}

Hope it helps

[MiguelSanches]

May I ask, where did you get that equation from?
I don’t think that is correct…

If you deduced it yourself from the equation below I can assure you that the one you posted is not correct, You can even check with a simple or double pendulum equation

For a simple pendulum:

Your denominator does not depend on the mass

For a double pendulum:

Your denominator must include cosine functions, which are not present in the equation you posted

[Rafsing]

I derived this equation from the equation number 38:

[MiguelSanches]

I don’t think that is correct.

The last term also depends on the acceleration. Which if I am not mistaken takes into account the cosine function that is missing you. Notice that for a simple pendulum this last term is zero because j = k, hence phi is 0. But for a double pendulum, when j = 1, and k = 2, then the phi function became 1 and introduces the cosine functions.

[Rafsing]

How would you transform it then?