//===========================ファイル名Quark2==========
//４次のルンゲ-クッタ法を用いたQuark運動？シミュレーションVer.1.2
//　　　　　　　　作成：加藤徳善　1997.3.23
//====================================================

import java.applet.*;
import java.awt.*;
//import java.lang.*;
//import java.util.*;

public class Quark2 extends Applet implements Runnable
{
    Graphics g_this;
    Image offScreen;
    Graphics offGraphics;
    Dimension dim1;

    Thread th;
    
    Color col=new Color(107,100,92); 
    Color col2=new Color(137,130,122); 
    Color col3=new Color(250,250,0); 
    
       
  	//ルンゲ-クッタ関連配列  
	double x[][]=new double[10][3];
	double v[][]=new double[10][3];
	double dx[][][]=new double[10][10][3];
	double nrm[][][]=new double[10][10][3];
	double pow[][]=new double[10][10];
	double ik[][][]=new double[4][10][3];
	double im[][][]=new double[4][10][3];
	double r[][]=new double[10][10];
	double aa[][]=new double[10][10];
	double al[][]=new double[10][10];
	int ap[][]=new int[10][10];
	double b[]=new double[10];
	double c[]=new double[10];
	double d[][]=new double[10][3];
	double x0[][]=new double[10][3];
	double mass[]=new double[10];
	double x1[][]=new double[10][3];
	double v1[][]=new double[10][3];
	double v2[][]=new double[10][3];
	double vn[][]=new double[10][3];
	double w[][]=new double[10][3];
	double w1[][]=new double[10][3];
	double va[]=new double[10];
	double f[][]=new double[10][3];
	double coef[]={1.,2.,2.,1.};
	
	//力のスイッチ
	boolean swa=true;//距離に関係する力（万有引力，弾性力など）
	boolean swb=true;//速さに比例する力（粘性抵抗）
	boolean swc=true;//速さの二乗に比例する力（慣性抵抗）
	boolean swd=true;//一定の力（地上の重力など）

	int dim;//次元-1
	int n;//物体数-1
	double dt;//時間間隔
	double time;//時間

    int FONTSIZE1=12;
	int oldx[][]=new int[3][2];
	int oldnv[]=new int[4];

    Point pt0=new Point(220,200);
    int drawmode=2;//スレッドによる　1:描画　2:描画停止
    int drawinitflag=1;// 1:画面初期化  0:何もしない
    int stepflag=0;//1:ストロボモード　0:何もしない
	int drawvflag=1;//速度ベクトルを表示する
	int targetn=0; 
	int oldtargetn=0;   
    int DD=15;
    int dd=10;
    Color co[]={Color.red,Color.green,Color.blue};
    
    public void init(){
    	int i,k;
        Font f=new Font("TimesRoman",Font.PLAIN,FONTSIZE1);//フォント作
        setFont(f);//フォント選択
 
        setLayout(new BorderLayout());//ボーダーレイアウトに
        Panel p = new Panel();//パネル作製
        p.setLayout(new GridLayout(4,1));
        p.add(new Button("Start"));//部品をパネルに貼り付け
        p.add(new Button("Stop")); 
        p.add(new Button("Strobo"));
        p.add(new Button("Reset"));
        add("East", p);//パネルを右に貼り付け


        initrk();//ルンゲ-クッタ関係の初期化
        //calc();
        drawinitflag=1;
        resize(500,400);
        setBackground(col);
        dim1=this.size();//現在の表示画面の大きさを調べる
        if (offScreen == null) {
	        offScreen = createImage(dim1.width,dim1.height);
	        offGraphics = offScreen.getGraphics();
        }        
   }

	//加速度の計算
	public void acce(){
		int i,j,k;
		double fw,fw1;
	
		if(swa||swc){
			for(i=0;i<=n;i++){
				if(swa){
					for(j=0;j<=i;j++){
						fw=0.;
						for(k=0;k<=dim;k++){
							if(i!=j){
								dx[i][j][k]=x1[i][k]-x1[j][k];
								dx[j][i][k]=x1[j][k]-x1[i][k];
							}
							else
								dx[i][j][k]=x1[i][k]-x0[j][k];
							fw+=dx[i][j][k]*dx[i][j][k];
						}
						r[i][j]=Math.sqrt(fw);
					
						for(k=0;k<=dim;k++){
							if(r[i][j]!=0)
								nrm[i][j][k]=dx[i][j][k]/r[i][j];
							else
								nrm[i][j][0]=1.;
							if(i!=j)
								nrm[j][i][k]=-nrm[i][j][k];
						}
						r[i][j]-=al[i][j];
						r[j][i]=r[i][j];
						if(ap[i][j]==0)
							pow[i][j]=1.;
						if(ap[i][j]==-2)
							pow[i][j]=1./r[i][j]/r[i][j];
						if(ap[i][j]==1)
							pow[i][j]=r[i][j];
						if(ap[i][j]==2){
							if(r[i][j]>al[i][j]){
								pow[i][j]=r[i][j]-al[i][j];
							}
							else{
								pow[i][j]=0.;
							}
						}
						pow[i][j]=pow[i][j]*aa[i][j];
						pow[j][i]=pow[i][j];
					}
				}
				fw1=0;
				for(k=0;k<=dim;k++)
					fw1+=v1[i][k]*v1[i][k];
				va[i]=Math.sqrt(fw1);
				for(k=0;k<=dim;k++){
					if(va[i]!=0.)
						vn[i][k]=v1[i][k]/va[i];
				}
			}
		}			

		for(i=0;i<=n;i++){
			for(k=0;k<=dim;k++)
				f[i][k]=0.;//clear
		}
		for(i=0;i<=n;i++){
			for(k=0;k<=dim;k++){
				if(swa){
					for(j=0;j<=n;j++)
						f[i][k]+=pow[i][j]*nrm[i][j][k];
				}
				if(swb)
					f[i][k]+=b[i]*v1[i][k];
				if(swc)
					f[i][k]+=c[i]*va[i]*va[i]*vn[i][k];
				if(swd)
					f[i][k]+=d[i][k];
				f[i][k]/=mass[i];
			}
		}		
	}
	
	//ルンゲ-クッタ法
	public void rk(){
		int i,j,k,l;
		for(i=0;i<=n;i++){
			for(k=0;k<=dim;k++){
				x1[i][k]=x[i][k];
				v1[i][k]=v[i][k];
				w[i][k]=0.;
				w1[i][k]=0.;
			}
		}
		acce();
		for(l=0;l<=2;l++){
			for(i=0;i<=n;i++){
				for(k=0;k<=dim;k++){
					ik[l][i][k]=dt*v1[i][k];
					im[l][i][k]=dt*f[i][k];
					x1[i][k]=x[i][k]+ik[l][i][k]/coef[l+1];
					v1[i][k]=v[i][k]+im[l][i][k]/coef[l+1];
					w[i][k]+=ik[l][i][k]*coef[l];
					w1[i][k]+=im[l][i][k]*coef[l];
				}
			}
			acce();
		}
		for(i=0;i<=n;i++){
			for(k=0;k<=dim;k++){
				ik[l][i][k]=dt*v1[i][k];
				im[l][i][k]=dt*f[i][k];
				x1[i][k]=x[i][k]+(w[i][k]+ik[3][i][k])/6.;
				v1[i][k]=v[i][k]+(w1[i][k]+im[3][i][k])/6.;
			}
		}
		for(i=0;i<=n;i++){
			for(k=0;k<=dim;k++){
				x[i][k]=x1[i][k];
				v[i][k]=v1[i][k];
			}
		}
	}
	
	public void initrk(){
		int i,j;
		dim=1;//次元-1
		n=2;//物体数-1
		dt=.4;//計算での時間間隔

		//物体の初めの位置x[i][n]　iは物体番号,nは座標軸番号
		//   x                y                 z
		x[0][0]=10 ; x[0][1]=0. ;
		x[1][0]=-10 ; x[1][1]=0. ; 
		x[2][0]=0. ; x[2][1]=-17.3; 

		//固定点の位置（物体ごとに力を受ける固定点を設定できる）
		//   x0              y0                 z0
		x0[0][0]=0. ; x0[0][1]=0. ; 
		x0[1][0]=0. ; x0[1][1]=0. ; 
		x0[2][0]=0. ; x0[2][1]=0. ; 

		//物体の初速度v[i][n]　iは物体番号,nは座標軸番号	
		//  v x             v y                 vz
		v[0][0]=0. ; v[0][1]=0 ; 
		v[1][0]=0. ; v[1][1]=0. ; v[1][2]=0.;	
		v[2][0]=0. ; v[2][1]=0. ; v[2][2]=0.;
	
		//質量
		mass[0]=2;
		mass[1]=2.;
		mass[2]=2.;
	
	//---力の設定-------	
		
		swa=true;//２物体間の相互作用
		//強さ
		aa[0][0]=-3.;
		aa[1][0]=-1.; aa[1][1]=-3.;
		aa[2][0]=-1.; aa[2][1]=-1.; aa[2][2]=-3.;
		//力のべき（-2は万有引力やクーロン力，1はバネの弾性力,2はゴムひもの弾性力）
		ap[0][0]=2;
		ap[1][0]=1; ap[1][1]=2;
		ap[2][0]=1; ap[2][1]=1; ap[2][2]=2;
		//バネの長さ
		al[0][0]=100.;
		al[1][0]=0.; al[1][1]=100.;
		al[2][0]=0.; al[2][1]=0.; al[2][2]=100.;
	
		if(n>0){//対称化（作用・反作用の法則）
			for(i=1;i<=n;i++){
				for(j=0;j<i;j++){
					aa[j][i]=aa[i][j];
					ap[j][i]=ap[i][j];
					al[j][i]=al[i][j];
				}
			}
		}
	
		swb=false;//速さに比例する力（粘性抵抗）
		//b[0]=0.; b[1]=0.; b[2]=0.;
		
		swc=false;//速さの二乗に比例する抵抗（慣性抵抗）
		//c[0]=0.; c[1]=0.; c[2]=-0.001.;	
		
		swd=false;//一定の力（重力など）
		//d[0][0]=0.; d[0][1]=9.8*mass[0];
		//d[1][0]=0.; d[1][1]=9.8*mass[1];
		//d[2][0]=0.; d[2][1]=9.8*mass[2];		
	}

     public void calc(){
 		rk();//ルンゲ-クッタ呼び出し

   }

		public boolean mouseDrag(Event e,int xd,int yd){
			int i;
			if(!((yd>380)&&(xd>360)&&(xd<440))){
			i=targetn;
			repaint();
			v[i][0]=(double)(xd-x[i][0]-pt0.x)/2.;
			v[i][1]=(double)(yd-x[i][1]-pt0.y)/2.;
			}					
			return true;
		}
		
		public boolean mouseDown(Event e,int xd,int yd){
			int i;
			if((yd>380)&&(xd>365)&&(xd<435)){
				oldtargetn=targetn;
				targetn=(xd-370)/20;
				repaint();
			}				
			return true;
		}
				
   
        public boolean action(Event e, Object o) {
        if(e.target instanceof Button) {
             if ("Start".equals(o)) {
                drawinitflag=1;
                stepflag=0;
                drawmode = 2;
                calc();
                drawmode = 1;
             } 
             if ("Stop".equals(o)) {
                 stepflag=0;
                 drawmode = 2;
                 return true;
             }
             if ("Strobo".equals(o)) {
                 drawinitflag=1;
                 stepflag=1;
                 drawmode = 1;
                 return true;
             }
             if ("Reset".equals(o)) {
                 drawinitflag=1;
                 stepflag=0;
                 drawmode = 2;
                 initrk();
                 drawvflag=1;
                 repaint();
                 return true;
             }
             return true;
      }
      return super.handleEvent(e);
   }

    public void start(){
        if(th==null){
             th=new Thread(this);
             th.start();
        }
    }
    public void run(){
        while(th!=null){
             if(drawmode==1){
                 try{
                      calc();
                      repaint();
                      th.sleep(50);
                 }catch(InterruptedException e){
                       ;
                 }
              }
        }
    }
    public void stop(){
        if(th !=null){
            th.stop();
            th=null;
        }
    }
    public void update(Graphics g){
        paint(g);//updateを無効にし、paintで置き換え
    }
    public void paint(Graphics g){
    	int i,j,k,d;
        if(drawinitflag==1){
             offGraphics.setColor(getBackground());
             offGraphics.fillRect(0,0,dim1.width,dim1.height);
             offGraphics.setColor(col2);
             offGraphics.fillOval(pt0.x-200,pt0.y-200,400,400);
           //背景色で画面全体を塗りつぶす（初期化）
            drawinitflag=0;
        }
         if((stepflag !=1)||(drawmode==2)){
	       	  offGraphics.setColor(getBackground());
             offGraphics.fillRect(0,0,dim1.width,dim1.height);
             offGraphics.setColor(col2);
             offGraphics.fillOval(pt0.x-200,pt0.y-200,400,400);
			}
			
		d=(int)(dd*.5);
		
        for(i=0;i<=n;i++){
        	oldx[i][0]=pt0.x+(int)(x[i][0]-dd*.5);
      		oldx[i][1]=pt0.y+(int)(x[i][1]-dd*.5);
      	}
      		
      	offGraphics.setColor(col3);
      	offGraphics.drawLine(oldx[0][0]+d,oldx[0][1]+d,oldx[1][0]+d,oldx[1][1]+d);
      	offGraphics.drawLine(oldx[2][0]+d,oldx[2][1]+d,oldx[1][0]+d,oldx[1][1]+d);
      	offGraphics.drawLine(oldx[0][0]+d,oldx[0][1]+d,oldx[2][0]+d,oldx[2][1]+d);
	
      	for(i=0;i<=n;i++){	
      		offGraphics.setColor(co[i]);
      		offGraphics.fillOval(oldx[i][0],oldx[i][1],(int)dd,(int)dd);
      		//if(drawinitflag==1)
      		offGraphics.fillOval(370+20*i,385,dd+1,dd+1);	
      		
      	}
      	
      	   
      		
      	//if((drawinitflag==1)||(oldtargetn!=targetn)){
//	      	offGraphics.setColor(co[oldtargetn]);
//			offGraphics.fillOval(370+20*oldtargetn,388,dd-5,dd-5);
			offGraphics.setColor(getBackground());
			offGraphics.fillOval(373+20*targetn,388,dd-5,dd-5);
	      	offGraphics.setColor(Color.white);
	      	offGraphics.drawString("Color selector",370,380);
// }
	      drawinitflag=0;
    	/*if(drawvflag==1){
			offGraphics.setColor(Color.pink);
  			oldnv[0]=pt0.x+(int)(x[0][0]);
  			oldnv[1]=pt0.y+(int)(x[0][1]);
  			oldnv[2]=pt0.x+(int)(x[0][0]+v[0][0]*2.);
  			oldnv[3]=pt0.y+(int)(x[0][1]+v[0][1]*2.);
  			offGraphics.drawLine(oldnv[0],oldnv[1],oldnv[2],oldnv[3]);
  		}*/
		g.drawImage(offScreen,0,0,this);
   }
}