// A circular ripple 
// Android Version of the Lattice Boltzmann Method
// Copyright @2013- Takeshi Seta All Rights Reserved.

package jp.seta.lbm_touch;

import android.os.Bundle;
import android.app.Activity;
import android.view.View;
import android.graphics.*;
import android.content.Context;
import android.os.Handler;
import android.os.Message;
import android.view.MotionEvent;

public class MainActivity extends Activity {
    LBMDrawView view;
    @Override
    public void onCreate(Bundle savedInstanceState) {
        super.onCreate(savedInstanceState);
        LBMDrawView v = new LBMDrawView(getApplication());
        setContentView(v);
    }
}

//Animation
class RedrawHandler extends Handler {
    private View view;
    private int delayTime;
    private int frameRate;
    public RedrawHandler(View view, int frameRate) {
        this.view = view;
        this.frameRate = frameRate;
    }
    public void start() {
        this.delayTime = 1000 / frameRate;
        this.sendMessageDelayed(obtainMessage(0), delayTime);
    }
    public void stop() {
        delayTime = 0;
    }
    @Override
    public void handleMessage(Message msg) {
        view.invalidate();
        if (delayTime == 0) return; // stop
        sendMessageDelayed(obtainMessage(0), delayTime);
    }
}

class LBMDrawView extends View {
    public int nx = 64, ny = 96, time = 0, in, jn, tx, ty, tflag = 0;
    public float tau = 0.8f, mu, tmp, u2, wid = 7f, u0 = 0.1f;
    float[][] u = new float[nx][ny], v = new float[nx][ny], rho = new float[nx][ny];
    float[][] x = new float[nx][ny], y = new float[nx][ny];
    float[][][] f = new float[9][nx][ny], f0 = new float[9][nx][ny], ftmp = new float[9][nx][ny];
    float[] cx = new float[9], cy = new float[9];
    Paint p = new Paint();
    Paint t = new Paint();

    public boolean onTouchEvent(MotionEvent event) {
        tx = (int)event.getX();
        ty = (int)event.getY();
        tx = (int)(tx/wid);
        ty = (int)(ty/wid);
        tflag = 1;
        return true;
    }    

    public LBMDrawView(Context c) {
        super(c);
        // initial condition
        mu = (2f*tau - 1f)/6f;
        for (int i = 0; i < nx; i++) { for (int j = 0; j < ny; j++) {
            u[i][j] = 0.0f;	v[i][j] = 0.0f; rho[i][j] = 1.0f;
            x[i][j] = i*wid; y[i][j] = j*wid;
        } }

        cx[0] = 0f; cy[0] = 0f;
        cx[1] = 1f; cy[1] = 0f; cx[2] = 0f; cy[2] = 1f;
        cx[3] =-1f; cy[3] = 0f; cx[4] = 0f; cy[4] =-1f;
        cx[5] = 1f; cy[5] = 1f; cx[6] =-1f; cy[6] = 1f;
        cx[7] =-1f; cy[7] =-1f; cx[8] = 1f; cy[8] =-1f;

        for (int i = 0; i < nx; i++) { for (int j = 0; j < ny; j++) {
            u2 = u[i][j]*u[i][j] + v[i][j]*v[i][j];      
            f[0][i][j] = rho[i][j]*(1f - 3f/2f*u2)*4f/9f;
            for (int k = 1; k < 5; k++) {
                tmp = cx[k]*u[i][j] + cy[k]*v[i][j];      
                f[k][i][j] = rho[i][j]*(1f + 3f*tmp + 9f/2f*tmp*tmp - 3f/2f*u2)/9f;
            }
            for (int k = 5; k < 9; k++) {
                tmp = cx[k]*u[i][j] + cy[k]*v[i][j];      
                f[k][i][j] = rho[i][j]*(1f + 3f*tmp + 9f/2f*tmp*tmp - 3f/2f*u2)/36f;
            }
         } }    

        setFocusable(true);
        RedrawHandler handler = new RedrawHandler(this, 8);
        handler.start();        
    }

    
    public void onDraw(Canvas c) { 
	//collision
        for (int i = 0; i < nx; i++) { for (int j = 0; j < ny; j++) {
            u2 = u[i][j]*u[i][j] + v[i][j]*v[i][j];      
            f0[0][i][j] = rho[i][j]*(1f - 3f/2f*u2)*4f/9f;
            for (int k = 1; k < 5; k++) {
                tmp = cx[k]*u[i][j] + cy[k]*v[i][j];      
                f0[k][i][j] = rho[i][j]*(1f + 3f*tmp + 9f/2f*tmp*tmp - 3f/2f*u2)/9f;
            }
            for (int k = 5; k < 9; k++) {
                tmp = cx[k]*u[i][j] + cy[k]*v[i][j];      
                f0[k][i][j] = rho[i][j]*(1f + 3f*tmp + 9f/2f*tmp*tmp - 3f/2f*u2)/36f;
            }
        } }
        for (int i = 0; i < nx; i++) { for (int j = 0; j < ny; j++) { for(int k = 0; k < 9; k++){
            f[k][i][j] = f[k][i][j] - (f[k][i][j] - f0[k][i][j])/tau;
        } } }

        // streaming
        for (int i = 0; i < nx; i++){ for (int j = 0; j < ny; j++){ for (int k = 0; k < 9; k++){
            ftmp[k][i][j] = f[k][i][j];
        } } }
        for (int i = 0; i < nx-1; i++){ for (int j = 0; j < ny; j++){
            in = i + 1;
            f[1][in][j] = ftmp[1][i][j];
        } }
        for (int i = 0; i < nx; i++){ for (int j = 0; j < ny-1; j++){
            jn = j + 1;
            f[2][i][jn] = ftmp[2][i][j];
        } }
        for (int i = 1; i < nx; i++){ for (int j = 0; j < ny; j++){
            in = i - 1;
            f[3][in][j] = ftmp[3][i][j];
        } }
        for (int i = 0; i < nx; i++){ for (int j = 1; j < ny; j++){
            jn = j - 1;
            f[4][i][jn] = ftmp[4][i][j];
        } }
        for (int i = 0; i < nx-1; i++){ for (int j = 0; j < ny-1; j++){
            in = i + 1;  jn = j + 1;
            f[5][in][jn] = ftmp[5][i][j];
        } }
        for (int i = 1; i < nx; i++){ for (int j = 0; j < ny-1; j++){
            in = i - 1; jn = j + 1;
            f[6][in][jn] = ftmp[6][i][j];
        } }
        for (int i = 1; i < nx; i++){ for (int j = 1; j < ny; j++){
            in = i - 1; jn = j - 1;
            f[7][in][jn] = ftmp[7][i][j];
        } }
        for (int i = 0; i < nx-1; i++){ for (int j = 1; j < ny; j++){
            in = i + 1;  jn = j - 1;
            f[8][in][jn] = ftmp[8][i][j];
        } }

        //boundary condition
        for (int j = 0; j < ny; j++){
            f[1][   0][j] = f[3][   0][j];
            f[5][   0][j] = f[7][   0][j];
            f[8][   0][j] = f[6][   0][j];
            f[3][nx-1][j] = f[1][nx-1][j];
            f[7][nx-1][j] = f[5][nx-1][j];
            f[6][nx-1][j] = f[8][nx-1][j];
         }
        for (int i = 0; i < nx; i++){
            f[2][i][   0] = f[4][i][   0];
            f[5][i][   0] = f[7][i][   0];
            f[6][i][   0] = f[8][i][   0];
            f[4][i][ny-1] = f[2][i][ny-1];
            f[7][i][ny-1] = f[5][i][ny-1];
            f[8][i][ny-1] = f[6][i][ny-1];
        }

        //physics
        for (int i = 0; i < nx; i++) { for (int j = 0; j < ny; j++) {
            rho[i][j] = f[0][i][j]; u[i][j] = 0f; v[i][j] = 0f;
            for (int k = 1; k < 9; k++){
                rho[i][j] = rho[i][j] + f[k][i][j];
                  u[i][j] =   u[i][j] + f[k][i][j] * cx[k];
                  v[i][j] =   v[i][j] + f[k][i][j] * cy[k];
            }
            if(rho[i][j] != 0){
                u[i][j] = u[i][j]/rho[i][j];
                v[i][j] = v[i][j]/rho[i][j];
            }else{
                u[i][j] = 0f; v[i][j] = 0f;
            }
        } }

        // A circular ripple 
        if(tflag == 1){
            for (int i = 0; i < nx; i++) { for (int j = 0; j < ny; j++) {
                if(i >=0){ if(i <nx){ if(j >=0){ if(j <ny){ if((i-tx)*(i-tx)+(j-ty)*(j-ty) !=0){
                //solitary wave
                u[i][j] = ((float)Math.exp(-(((i-tx)*(i-tx)+(j-ty)*(j-ty))/(4f*1f)))/(2f*(float)Math.sqrt(Math.PI*1f))*(float)(i-tx)/((float)Math.sqrt((i-tx)*(i-tx)+(j-ty)*(j-ty))));
                v[i][j] = ((float)Math.exp(-(((i-tx)*(i-tx)+(j-ty)*(j-ty))/(4f*1f)))/(2f*(float)Math.sqrt(Math.PI*1f))*(float)(j-ty)/((float)Math.sqrt((i-tx)*(i-tx)+(j-ty)*(j-ty))));
                //soliton
//              u[i][j] = ((2f*0.5f*0.5f)*(float)(i-tx)/((float)Math.sqrt((i-tx)*(i-tx)+(j-ty)*(j-ty))))/((float)Math.cosh(0.5f*(Math.sqrt((i-tx)*(i-tx)+(j-ty)*(j-ty)+0f)))*(float)Math.cosh(0.5f*(Math.sqrt((i-tx)*(i-tx)+(j-ty)*(j-ty)+0f))));
//              v[i][j] = ((2f*0.5f*0.5f)*(float)(j-ty)/((float)Math.sqrt((i-tx)*(i-tx)+(j-ty)*(j-ty))))/((float)Math.cosh(0.5f*(Math.sqrt((i-tx)*(i-tx)+(j-ty)*(j-ty)+0f)))*(float)Math.cosh(0.5f*(Math.sqrt((i-tx)*(i-tx)+(j-ty)*(j-ty)+0f))));  
                } } } } }
            } }
            tflag = 0;
        }
        tmp = 0f;
        for (int i = 0; i < nx; i++) { for (int j = 0; j < ny; j++) {
            u2 = (float)Math.sqrt(u[i][j]*u[i][j] + v[i][j]*v[i][j]);
            if(u2 > tmp){ tmp = u2; }
        } }


        //graphics
        c.drawColor(Color.WHITE);
        p.setStrokeWidth(3);
        p.setStyle(Paint.Style.STROKE);
        p.setColor(Color.GRAY);
        c.drawRect(x[0][0], y[0][0], x[nx-1][ny-1]+ wid, y[nx-1][ny-1]+ wid, p);            
        p.setStrokeWidth(2);
        for (int i = 0; i < nx; i++) { for (int j = 0; j < ny; j++) {
            u2 = (float)Math.sqrt(u[i][j]*u[i][j] + v[i][j]*v[i][j]);      
            p.setColor(Color.argb((int)(u2/u0*255),0,0,(int)(u2/u0*255)));
            p.setStyle(Paint.Style.FILL);
            c.drawRect(x[i][j], y[i][j], x[i][j] + wid, y[i][j] + wid, p); 
            p.setColor(Color.argb(255,0, 255, 255));
            c.drawLine(x[i][j] + wid/2, y[i][j] + wid/2, x[i][j] + wid/2 + u[i][j]/u0*wid*10, y[i][j] + wid/2 + v[i][j]/u0*wid*10, p); 
        } }

        // Text
        t.setColor(Color.RED);
        t.setTextSize (16);
        String msg1 = "Tap the screen";
        c.drawText(msg1, 10, 30, t);

        t.setColor(Color.BLACK);
        String msg2 = "time: " +time++ +" , ( " + tx +" ," +ty+" ) " +", velocity " + tmp;
        c.drawText(msg2, 130, 30, t);
        String msg3 = "Android Version of the Lattice Boltzmann Method";
        c.drawText(msg3, 10, 700, t);
        String msg4 = "Copyright @2013- Takeshi Seta All Rights Reserved.";
        c.drawText(msg4, 10, 720, t);        
    }
}