// Program name - harmony.java - PhD Project
// written by Jeff Mirick
// February 14, 1999
//  There are 3 classes in the program
//  1.)  myfunction - defines the evolution of the
//       harmonic oscillator in time.
//  2.)  plot - Defines the ground sate for the harmonic oscillator problem
//       for the two values of k.
//  3.)  harmony - Is the main executable class defining the animation,
//       panels, and graphical interfaces.

import java.awt.*;
import java.awt.event.*;
import java.util.*;
import java.applet.*;
import java.awt.image.ImageObserver;

class myfunction {
   int i;
   double w_not,pi;
   double y[] = new double[500];
   double potential[] = new double[500];
   double fk1[] = new double[500];
   double xx;
   int maxx,maxy,minx,miny;
   double xlow,xhigh,ylow,yhigh,E0;
   public void init() {
   	xlow=-3.0;
      xhigh=3.0;
      ylow=0.0;
      yhigh=3.0;
      minx=0;
      miny=0;
      maxx=350;
		maxy=300;
      pi=3.14159265;
    }

	void compute(int itime, double k1, double k2) {
      //
   	//  Defines time evolution of harmonic oscillator
   	//
      double a127_2,a125_2,a123_2,a12_2,cos2,cos4,cos6;
      double a1,a2,a12,a123,a125,c1,e,prob;
      double psi0,psi2,psi4,psi6,sa1a2,a2p,t,factor;
      w_not=Math.sqrt(k2);
   	a1=Math.sqrt(k1);
   	a2=Math.sqrt(k2);
   	sa1a2=Math.sqrt(Math.sqrt(a1*a2));
   	a12=a1+a2;
   	a2p=Math.sqrt(Math.sqrt(a2/pi));
   	a123=a12*a12*a12;
   	a125=a123*a12*a12;
   	c1=a2p*sa1a2;
      a12_2=Math.sqrt(2.0/a12);
      a123_2=Math.sqrt(2.0/a123);
      a125_2=Math.sqrt(2.0/a125);
      a127_2=Math.sqrt(2.0/a125/a12/a12);
      //
   	//
   	t=(double)(itime)/50.0;
      cos2=Math.cos(2.0*w_not*t);
      cos4=Math.cos(4.0*w_not*t);
      cos6=Math.cos(6.0*w_not*t);
     	for(i=0;i<maxx;i++) {
         xx=xlow+(double)(i)*(xhigh-xlow)/(maxx-minx);
         e=Math.exp(-a2*xx*xx/2.0);
      	psi0=c1*a12_2*e;
      	psi2=c1*(2.0*a2*a123_2-a12_2)*(2.0*a2*xx*xx-1.0)*e/2.0;
    		psi4=c1*(6.0*a2*a2*a125_2-6.0*a2*a123_2+
         		1.5*a12_2)*(2.0*a2*a2*xx*xx*xx*xx
               	-6.0*a2*xx*xx+1.5)*e/6.0;
         psi6=c1*(60.0*a2*a2*a2*a127_2-90.0*a2*a2*a125_2+45.0*a2*a123_2
                -7.5*a12_2)*(8.0*a2*a2*a2*xx*xx*xx*xx*xx*xx
                -60.0*a2*a2*xx*xx*xx*xx+90.0*a2*xx*xx-15.0)*e/360.0;
      	prob=psi0*psi0+psi2*psi2+psi4*psi4+
      			2.0*(psi0*psi2*cos2+psi0*psi4*cos4+psi2*psi4*cos2);
         prob=prob+psi6*psi6+2.0*(psi0*psi6*cos6+psi2*psi6*cos4+
                         psi4*psi6*cos2);
     		y[i]=prob+Math.sqrt(k1)/2.0;
      }
	}

   void f3(double k1)
	{
    for(i=0;i<maxx;i++) {
         xx=xlow+(double)(i)*(xhigh-xlow)/(maxx-minx);
    		potential[i]=k1*xx*xx/2.0;
    }
	}

     void f1(double k1) {
   	double a1p,a1,x1;
      E0=Math.sqrt(k1)/2.0;
   	a1=Math.sqrt(k1);
      a1p=Math.sqrt(a1/pi);
      for(i=0;i<maxx;i++) {
      	xx=xlow+(double)(i)*(xhigh-xlow)/(maxx-minx);
         x1=a1*xx*xx;
   		fk1[i]= a1p*Math.exp(-x1)+E0;
      }
	}

    public void paint(Graphics g, ImageObserver obs) {
    	int ix,iy,ix1,iy1;
      g.setColor(Color.blue);
      ix=minx;
   	iy=maxy-2-(int)((potential[0]-ylow)*(maxy-miny)/(yhigh-ylow));
   	for(i=1;i<maxx;i++) {
      	xx=xlow+(double)(i)*(xhigh-xlow)/(maxx-minx);
			ix1=minx+(int)((xx-xlow)*(double)(maxx-minx)/(xhigh-xlow));
			iy1=maxy-2-(int)((potential[i]-ylow)*(maxy-miny)/(yhigh-ylow));
      	if (iy1>miny) g.drawLine(ix,iy,ix1,iy1);
      	iy=iy1;	ix=ix1;
		}

      g.setColor(Color.green);
       g.drawLine(0, maxy, maxx, maxy);
      ix=minx;
   	iy=maxy-2-(int)((y[0]-ylow)*(maxy-miny)/(yhigh-ylow));
   	for(i=1;i<maxx;i++) {
      	xx=xlow+(double)(i)*(xhigh-xlow)/(maxx-minx);
			ix1=minx+(int)((xx-xlow)*(double)(maxx-minx)/(xhigh-xlow));
			iy1=maxy-2-(int)((y[i]-ylow)*(maxy-miny)/(yhigh-ylow));
      	g.drawLine(ix,iy,ix1,iy1);
      	iy=iy1;	ix=ix1;
		}

      g.setColor(Color.yellow);
       g.drawLine(0, maxy, maxx, maxy);
      ix=minx;
   	iy=maxy-2-(int)((fk1[0]-ylow)*(maxy-miny)/(yhigh-ylow));
   	for(i=1;i<maxx;i++) {
      	xx=xlow+(double)(i)*(xhigh-xlow)/(maxx-minx);
			ix1=minx+(int)((xx-xlow)*(double)(maxx-minx)/(xhigh-xlow));
			iy1=maxy-2-(int)((fk1[i]-ylow)*(maxy-miny)/(yhigh-ylow));
      	g.drawLine(ix,iy,ix1,iy1);
      	iy=iy1;	ix=ix1;
		}
      g.setColor(Color.red);
      iy=maxy-2-(int)((E0-ylow)*(maxy-miny)/(yhigh-ylow));
      g.drawLine(minx,iy,maxx,iy);
 	}
}

class plot {
   int i,minx,miny,maxx,maxy;
   double k1,k2,x1,x2,xx;
   double pi,w_not;
   double xlow,xhigh,ylow,yhigh;
   double fk1[] = new double[500];
   double fk2[] = new double[500];
   double potential[] = new double[500];
   public void init() {
      xlow=-3.0;
      xhigh=3.0;
      ylow=0.0;
      yhigh=3.0;
      minx=0;
      miny=0;
      maxx=350;
		maxy=300;
      pi=3.14159265;
      }

      void f1(double k1) {
      double a1p,a1,x1;
      w_not=Math.sqrt(k1);
     	a1=Math.sqrt(k1);
   	a1p=Math.sqrt(a1/pi);
      for(i=0;i<maxx;i++) {
      	xx=xlow+(double)(i)*(xhigh-xlow)/(maxx-minx);
         x1=a1*xx*xx;
      	fk1[i]=a1p*Math.exp(-x1)+w_not/2.0;
    	}
    }

   void f2(double k2) {
   	double a2p,a2,x2;
   	a2=Math.sqrt(k2);
      a2p=Math.sqrt(a2/pi);
      for(i=0;i<maxx;i++) {
      	xx=xlow+(double)(i)*(xhigh-xlow)/(maxx-minx);
         x2=a2*xx*xx;
   		fk2[i]= a2p*Math.exp(-x2)+w_not/2.0;
      }
	}

   void f3(double k1)
	{
   	for(i=0;i<maxx;i++) {
      	xx=xlow+(double)(i)*(xhigh-xlow)/(maxx-minx);
   		potential[i]=k1*xx*xx/2.0;
      }
	}
	public void paint(Graphics g, ImageObserver obs) {

   	int ix,iy,ix1,iy1;

   	g.setColor(Color.yellow);
      g.drawLine(0, maxy, maxx, maxy);

      // Plot of function for k1

      ix=minx;
   	iy=maxy-2-(int)((fk1[0]-ylow)*(maxy-miny)/(yhigh-ylow));
   	for(i=1;i<maxx;i++) {
      	xx=xlow+(double)(i)*(xhigh-xlow)/(maxx-minx);
			ix1=minx+(int)((xx-xlow)*(double)(maxx-minx)/(xhigh-xlow));
			iy1=maxy-2-(int)((fk1[i]-ylow)*(maxy-miny)/(yhigh-ylow));
      	g.drawLine(ix,iy,ix1,iy1);
      	iy=iy1;	ix=ix1;
		}

      // Plot of function for k2

      g.setColor(Color.green);
      ix=minx;
   	iy=maxy-2-(int)((fk2[0]-ylow)*(maxy-miny)/(yhigh-ylow));
   	for(i=1;i<maxx;i++) {
      	xx=xlow+(double)(i)*(xhigh-xlow)/(maxx-minx);
			ix1=minx+(int)((xx-xlow)*(double)(maxx-minx)/(xhigh-xlow));
			iy1=maxy-2-(int)((fk2[i]-ylow)*(maxy-miny)/(yhigh-ylow));
      	if (iy1>miny) g.drawLine(ix,iy,ix1,iy1);
      	iy=iy1;	ix=ix1;
		}

      g.setColor(Color.blue);
      ix=minx;
   	iy=maxy-2-(int)((potential[0]-ylow)*(maxy-miny)/(yhigh-ylow));
   	for(i=1;i<maxx;i++) {
      	xx=xlow+(double)(i)*(xhigh-xlow)/(maxx-minx);
			ix1=minx+(int)((xx-xlow)*(double)(maxx-minx)/(xhigh-xlow));
			iy1=maxy-2-(int)((potential[i]-ylow)*(maxy-miny)/(yhigh-ylow));
      	g.drawLine(ix,iy,ix1,iy1);
      	iy=iy1;	ix=ix1;
		}

      g.setColor(Color.red);
      iy=maxy-2-(int)((w_not/2.0-ylow)*(maxy-miny)/(yhigh-ylow));
      g.drawLine(minx,iy,maxx,iy);
   }
}


public class harmony extends Applet
		implements ItemListener, ActionListener, Runnable {
   int itime;
   double k1,k2;
   Panel mypanel1;
   Button but1, but2;
   Choice mychoice1, mychoice2;
	Image offimg;
   Dimension offsize;
   Graphics offg;
   Thread myanimation;
	myfunction myharmony = new myfunction();
   plot myplot = new plot();

   public void init() {
   	setBackground(Color.black);
      setForeground(Color.green);
     	setLayout(new BorderLayout());

   	mypanel1 = new Panel();
  		mychoice1=new Choice();
  		mychoice2=new Choice();

      but1=new Button("Start Program");
      but1.addActionListener(this);
      mypanel1.add(but1);

     	mychoice1.addItem("1.0");
 		mychoice1.addItem("2.0");
 		mychoice1.addItem("3.0");
      mychoice1.addItem("4.0");
      mychoice1.addItem("5.0");
 		mychoice1.addItem("6.0");
      mychoice1.addItem("7.0");
      mychoice1.addItem("8.0");
      mychoice1.addItem("9.0");
      mychoice1.addItem("10.0");
 		mypanel1.add(mychoice1);

 		mychoice2.addItem("1.0");
 		mychoice2.addItem("2.0");
 		mychoice2.addItem("3.0");
 		mychoice2.addItem("4.0");
      mychoice2.addItem("5.0");
 		mychoice2.addItem("6.0");
      mychoice2.addItem("7.0");
      mychoice2.addItem("8.0");
      mychoice2.addItem("9.0");
      mychoice2.addItem("10.0");
		mypanel1.add(mychoice2);

 		add("South",mypanel1);

      mychoice1.addItemListener(this);
      mychoice2.addItemListener(this);
   	itime=0;
   	k1=1.0;
   	k2=1.0;
    	myharmony.init();
    	myplot.init();
    	myplot.f1(k1);
    	myplot.f2(k2);
    	myplot.f3(k1);
    	offsize=getSize();
    	offimg=createImage (offsize.width, offsize.height);
    	offg=offimg.getGraphics();
   }

   public void stop(){
   	if(myanimation!=null) {
	 	myanimation=null;
      but1.setLabel("Start Animation");
      }
      }

	public void run() {
   	Thread thisThread=Thread.currentThread();
   	offsize=getSize();
      offimg=createImage (offsize.width, offsize.height);
      offg=offimg.getGraphics();
      itime=0;
      while (myanimation==thisThread) {
    		itime=itime+1;
         myharmony.f3(k1);
         myharmony.f1(k1);
		 	myharmony.compute(itime,k1,k2);
   		repaint();
			try { Thread.sleep(100); }
   		catch ( InterruptedException e) {break; }
   	}
	}

   public synchronized void update (Graphics g) {
	  	offg.setColor(Color.black);
      offg.fillRect(0,0,offsize.width, offsize.height);
      if(but1.getLabel()=="Start Animation") myplot.paint(offg, this);
      else myharmony.paint(offg,this);
      g.drawImage(offimg, 0, 0, this);
   }

   public void actionPerformed(ActionEvent event) {
   	Object object = but1.getLabel();
      if (object=="Start Program") {
      	but1.setLabel("Start Animation");
         repaint();
      }
      if (object=="Start Animation") {
      	but1.setLabel("Stop Animation");
         myanimation = new Thread(this);
         myanimation.start();
      }
      else {
        	but1.setLabel("Start Animation");
         myanimation=null;
         repaint();
       }
   }
	public void itemStateChanged(ItemEvent event) {
   	String str1, str2;
      str1=mychoice1.getSelectedItem();
      str2=mychoice2.getSelectedItem();
      k1=Double.valueOf(str1).doubleValue();
      k2=Double.valueOf(str2).doubleValue();
      myplot.f1(k1);
      myplot.f2(k2);
      myplot.f3(k1);
      repaint();
   }
}