Physics 310 - Intermediate Mechanics


Solving problems numerically

You may need to write simple computer programs to solve some problems. Windows is not particularly well suited for doing this. For this reason, here are instructions for writing simple programs and compiling them on one of the Sun workstations that use UNIX for their operating system.

Logging in

You can find Sun workstations in the physics building and log in directly using your ICS account. Alternatively, you can log in to any Windows computer and connecto to a Sun using ssh. To do this you need to do this:
Start --> All Programs
      --> Standard Software
      --> Telecommunications
      --> Secure CRT4.1
      --> Secure CRT4.1
From here you can log in to expert.ics.purdue.edu using your ICS account.

Editing your program

You first need to type in the program you need to develop to solve the problem on which you are working. There are several editors available: Once you have logged in, if the backspace key doesn't work, try this:
expert.ics.purdue.edu% stty erase '^?'
Once you have logged in you can create a new directory in which to do your work, and edit a new program like this:
expert.ics.purdue.edu% mkdir phys310
expert.ics.purdue.edu% cd phys310
expert.ics.purdue.edu% pico junk.c
Then you can enter a program like this:
#include "stdio.h"

main() {
  printf( "Hello world.\n" );
}
Then you can compile and run it like this:
expert.ics.purdue.edu% cc junk.c
expert.ics.purdue.edu% ./a.out
Hello world.

A better example

Most of the problems for physics 310 involve mathematical operations like sin(), cos(), etc. These are contained in a math library. Here is a more realistic program, that solves the differential equation:
dy/dt = -y t
using the 4-th order Runge-Kutta method and compares it with the exact solution:

#include "stdio.h"
#include "math.h"

double dydt(double t,double y) {
  return -y*t;
}

main() {
  double y, y0, y1, yexact, h, t, tmax;
  double k1, k2, k3, k4;

  h = 0.1;
  y0 = y = 1;
  t = 0;
  tmax = 5;
  printf( "   t     y (numerical)     y (exact)\n" );
  printf( "--------  -------------  ------------\n" );
  while ( t < tmax ) {
    k1 = dydt(t,y);
    k2 = dydt(t+0.5*h,y+0.5*h*k1);
    k3 = dydt(t+0.5*h,y+0.5*h*k2);
    k4 = dydt(t+h,y+h*k3);
    y1 = y + h*(k1+2*k2+2*k3+k4)/6;
    y = y1;
    t = t + h;
    yexact = y0*exp(-0.5*t*t);
    printf( "%6f  %12f  %12f\n", t, y, yexact );
  }
}
You can compile and run this program like this:
expert.ics.purdue.edu% cc junk.c -lm
expert.ics.purdue.edu% ./a.out
   t     y (numerical)     y (exact)
-------  -------------  ------------
0.100000      0.995012      0.995012
0.200000      0.980199      0.980199
0.300000      0.955997      0.955997
0.400000      0.923116      0.923116
0.500000      0.882497      0.882497
0.600000      0.835270      0.835270
0.700000      0.782705      0.782705
0.800000      0.726149      0.726149
0.900000      0.666977      0.666977
1.000000      0.606531      0.606531
1.100000      0.546075      0.546074
1.200000      0.486752      0.486752
1.300000      0.429558      0.429557
1.400000      0.375311      0.375311
1.500000      0.324653      0.324652
1.600000      0.278038      0.278037
1.700000      0.235747      0.235746
1.800000      0.197900      0.197899
1.900000      0.164476      0.164474
2.000000      0.135337      0.135335
etc...

Other languages...

For a blast from the past, you might also enjoy programming in FORTRAN. Here is the equivalent source code:
      REAL FUNCTION DYDT(T,Y)
      IMPLICIT NONE
      REAL T, Y

      DYDT = -Y*T
      RETURN
      END

      PROGRAM RUNGE
      IMPLICIT NONE

      REAL Y, Y0, Y1, YEXACT, H, T, TMAX
      REAL K1, K2, K3, K4

      REAL DYDT

      H = 0.1
      Y0 = 1
      Y = 1
      TMAX = 5
      WRITE(6,'("   T       Y (NUMERICAL)      Y (EXACT)")')
      WRITE(6,'("--------   -------------   ------------")')
  100 CONTINUE
      K1 = DYDT(T,Y)
      K2 = DYDT(T+0.5*H,Y+0.5*H*K1)
      K3 = DYDT(T+0.5*H,Y+0.5*H*K2)
      K4 = DYDT(T+H,Y+H*K3)
      Y1 = Y + H*(K1+2*K2+2*K3+K4)/6
      Y = Y1
      T = T + H
      YEXACT = Y0*EXP(-0.5*T*T)
      WRITE(6,'(F8.6,2X,F12.6,2X,F12.6)') T, Y, YEXACT
      IF ( T .LT. TMAX ) GOTO 100
      STOP
      END
You can compile it and run it like this:
expert.ics.purdue.edu% f77 junk.f
expert.ics.purdue.edu% ./a.out
   T       Y (NUMERICAL)      Y (EXACT)
--------   -------------   ------------
0.100000      0.995012      0.995012
0.200000      0.980199      0.980199
0.300000      0.955997      0.955997
0.400000      0.923116      0.923116
0.500000      0.882497      0.882497
0.600000      0.835270      0.835270
0.700000      0.782705      0.782705
0.800000      0.726149      0.726149
etc...