Magnetic pendulum part 2

Jonathon Sumner
Physics Department
Dawson College
Jean-François Brière
Physics Department
Dawson College
Sameer Bhatnagar
Physics Department
Dawson College
March 27, 2019
March 27, 2019

From the magnetic pendulum to chaos

Now that we have a working magnetic pendulum model, let's try to find out where the magnet ends up based on its initial position.

Lab instructions

Magnetic pendulum instructions (complete)

Code for Task B

Magnetic pendulum map

Code sample: Magnetic pendulum map
/**
   Meeting 6 - Task B

   magneticPendulumMap

   Author: J.F. Briere, J. Sumner
   Current version written: March 2019
   Description: Map a square of the choatic behviour of the magnetic pendulum,
				returning a number that is assigned to the magnet where the pendulum
				is ending its trajectory
*/
import java.io.*;
import java.lang.*;
import java.util.Locale;

public class magneticPendulumMapStudents
{
	///////////////////////////////////////////////////////////////////////////////////////////////
	// 1. Declare constants
	// Remember to put units!
	///////////////////////////////////////////////////////////////////////////////////////////////
    public static final double DT = 1e-2;     // s
    public static final double Tmax = 60.;    // s

    // Parameters of the equations
    public static final double k = 0.5;       // spring constant N/m/kg
    public static final double b = 0.15;	  // drag coefficient Ns/m/kg
    public static final double z = 0.25;      // vertical distance between the pendulum and the magnets m
	public static final double J = 1.0;       // strenght of magnetic force in Nm^2/kg
 
	// Magnet location
    public static final double[] Mx = {1,-0.5, -0.5};	// X position of magnets m
    public static final double[] My = {0,Math.sqrt(3.)/2., -Math.sqrt(3.)/2.};	// Y position of magnets m

    // Set location and size of square within which we will map the behaviour of the chaotic pendulum
    // Note: You can use any value here from (-4.8,-4.8) to (4.,4)
    public static final double mapSquareLowerLeftCornerX = -0.; //Please use multiples of 0.8
    public static final double mapSquareLowerLeftCornerY = -0.; //Please use multiples of 0.8
    public static final double mapRangeX = 0.8; // Please do not change
    public static final double mapRangeY = 0.8; // Please do not change

	// Set grid size. Will calculate gridSize x gridSize points
	public static final int gridSize=40; //Size of grid in each direction

    // Start main method
    public static void main(String[] args)
    {
		//Stating a counter for time
		long startTime = System.currentTimeMillis();

		///////////////////////////////////////////////////////////////////////////////////
		// 2. Open file to write the data
		///////////////////////////////////////////////////////////////////////////////////
		PrintWriter outputFile = null;
		try
		{
		outputFile = new PrintWriter(new FileOutputStream("magnetic_pendulum_map.txt",false));
		}
		catch(FileNotFoundException e)
		{
		System.out.println("File error.  Program aborted.");
		System.exit(0);
		}

		///////////////////////////////////////////////////////////////////////////////////
		// 3. Creating the map
		///////////////////////////////////////////////////////////////////////////////////

		int magnet; // Numbers between 1 and 3 telling the final attracting magnet
		double xInitial, yInitial;
		double dx = mapRangeX/gridSize;
		double dy = mapRangeY/gridSize;
		
		for (int i=0; i<gridSize; i++)
		{
			xInitial = mapSquareLowerLeftCornerX + i*dx; //Initial x of pendulum
			for (int j=0; j<gridSize; j++)
			{
				yInitial = mapSquareLowerLeftCornerY + j*dy;//Initial y of pendulum

				magnet = trajectory(xInitial,yInitial);
				outputFile.printf(Locale.US,"%6.5f\t%6.5f\t%2d\n",xInitial,yInitial,magnet);
				//outputFile2.println(xInitial + "	    " + yInitial + "	" + magnet);
			}
			System.out.println("Done with line "+i+" of "+gridSize+". x = "+xInitial);
		}

		//Cleaning up
		outputFile.close();
		
		// Counting and printing running time.
		long endTime = System.currentTimeMillis();
		System.out.println("The code simulated "+ Math.pow(gridSize,2)+" trajectories in "+(endTime - startTime)/1000+" s");
    }

    ///////////////////////////////////////////////////////////////////////////////////
  	// 5. Calculating acceleration in x
	// parameters: 
	//	xpos = position of pendulum in x
	// 	ypos = position of pendulum in y-acceleration
	// 	xvel = x component of velocity of the pendulum
  	///////////////////////////////////////////////////////////////////////////////////
	public static double accelerationx(double xpos, double ypos, double xvel)
    {
		double accx;
		
		accx = -k*xpos -b*xvel;

		// TASK B.4: ADD YOUR CODE TO CALCULATE THE EFFECT
		// OF THE MAGNETS FROM TASK A

		return accx;
    }

    ///////////////////////////////////////////////////////////////////////////////////
  	// 6. Calculating acceleration in y
 	// parameters: 
	//	xpos = position of pendulum in x
	// 	ypos = position of pendulum in y-acceleration
	// 	yvel = x component of velocity of the pendulum
 	///////////////////////////////////////////////////////////////////////////////////
    public static double accelerationy(double xpos, double ypos, double yvel)
    {
		double accy;

		// Spring and damping contribution
		accy = -k*ypos -b*yvel;

		// TASK B.4: ADD YOUR CODE TO CALCULATE THE EFFECT
		// OF THE MAGNETS FROM TASK A
	
		return accy;
    }
	
	///////////////////////////////////////////////////////////////////////////////////
  	// 7. Calculating the pendulum trajectory
 	// parameters: 
	//	x0 = initial position of pendulum in x
	// 	y0 = initial position of pendulum in y
	// returns:	number of magnet where it ends up
 	///////////////////////////////////////////////////////////////////////////////////
	public static int trajectory(double x0, double y0)
    {
		// Calculate length of arrays
		int imax = (int)(Tmax/DT);		// Maximal index

		// Declare main variables
		double[] t = new double[imax];		// time in sec

		// Mass 1
		double[] x1 = new double[imax];  // x-position in m
		double[] vx1 = new double[imax];// x-velocity in m/s
		double[] ax1 = new double[imax];// x-acceleration in m/s/s

		double[] y1 = new double[imax];	// y-position in m
		double[] vy1 = new double[imax];// y-velocity in m/s
		double[] ay1 = new double[imax];// y-acceleration in m/s/s

		// Initialize the first point;
		t[0] = 0;

		x1[0] = x0;
		y1[0] = y0;

		vx1[0] = 0;
		vy1[0] = 0;

		ax1[0] = accelerationx(x1[0],y1[0],vx1[0]);
		ay1[0] = accelerationy(x1[0],y1[0],vy1[0]);

		// Do one step of Euler to get value for i=1
		t[1] = t[0]+DT;

		// Mass 1
		x1[1] = x1[0] + vx1[0]*DT;
		vx1[1] = (x1[1]-x1[0])/DT;
		ax1[1] = accelerationx(x1[1],y1[1],vx1[1]);

		y1[1] = y1[0] + vy1[0]*DT;
		vy1[1] = (y1[1]-y1[0])/DT;
		ay1[1] = accelerationy(x1[1],y1[1],vy1[1]);


		for(int i=2;i<imax;i++)
			{
			// Update time
			t[i] = t[i-1]+DT;

			// Update position
			x1[i] = 2*x1[i-1] - x1[i-2] + ax1[i-1]*Math.pow(DT,2);
			y1[i] = 2*y1[i-1] - y1[i-2] + ay1[i-1]*Math.pow(DT,2);

			// Update velocity
			vx1[i] = vx1[i-1] + ax1[i-1]*DT;
			vy1[i] = vy1[i-1] + ay1[i-1]*DT;

			// Update acceleration
			ax1[i] = accelerationx(x1[i],y1[i],vx1[i]);
			ay1[i] = accelerationy(x1[i],y1[i],vy1[i]);
		}

		// TASK B.5: Add a few lines that will determine near which magnet the mass is at
		// the end of the simulation.
		// Represent the magnets as integers from 1 to 3 and assign a 4 if the magnet is
		// not close to a magnet.
		int magnet = 0; 

		return magnet;
    }
}

To submit your square

Go to:

https://compute.dawsoncollege.qc.ca/magnetic_pendulum/

Codes for Task C

Two codes for Task C:

Codes to compare trajectories