//		-*- Mode: Java -*-
// Author(s):	J. Robert Buchanan and Zhoude Shao
// Address:	Department of Mathematics
//		Millersville University
//      	P.O. Box 1002
//		Millersville, PA 17551-0302
// Phone:	717-871-7305
// FAX:		717-871-7948
// Email:	Robert.Buchanan@millersville.edu, Zhoude.Shao@millersville.edu
//
// Date:	11/02/2016
//
// Purpose:	Explicit finite difference method for solving the
//		homogeneous wave equation with Dirichlet boundary conditions.
//		This program accompanies Exercise 14.4.2.
//
// To compile:	javac exercise14.java
// To execute:	java exercise14
//
import java.io.*;

public class exercise14
{
    public static double f(double x)	// initial displacement
    {
	double result	;
	
	result = 4.0*x*(1.0-x)	;

	return(result)	;
    }	// end of function f( )

    public static double g(double x)	// initial velocity
    {
	double result	;

	result = Math.sin(3.0*Math.PI*x)	;
	
	return(result)	;
    }	// end of function g( )

    public static void main(String args[])
    {
	int i, j;	// counter(s)
	int nelts = 10+1;
	double h;	// spatial discretization step size
	double k;	// temporal discretization step size
	double r;	// method parameter
	double [][] soln = new double[3][nelts];// results for j-1,j,j+1

	// Set step sizes
	h = 1.0/((double)(nelts-1))	;
	k = 0.05	;
	r = k/h	;
	// Set initial conditions
	for ( i=0 ; i<nelts ; i++ ) {
	    soln[0][i] = f(((double)i)*h)	;
	    if ( (i==0) || (i==nelts-1) ) soln[1][i] = 0.0	;
	    else {
		soln[1][i] = 0.5*r*r*f(((double)(i+1))*h)	;
		soln[1][i]+= (1.0-r*r)*f(((double)i)*h)	;
		soln[1][i]+= 0.5*r*r*f(((double)(i-1))*h)	;
		soln[1][i]+= k*g(((double)i)*h)	;
	    }	// end of else clause
	}	// end of for ( i ) loop
	// Output the first two time steps
	for ( i=0 ; i<nelts ; i++ ) System.out.printf("%8.4f ",soln[0][i]);
	System.out.printf("\n");
	for ( i=0 ; i<nelts ; i++ ) System.out.printf("%8.4f ",soln[1][i]);
	System.out.printf("\n");
	// Iterate for t = 0, k, 2k, ..., 10k
	for ( j=2 ; j<10+1 ; j++ ) {
	    for ( i=0 ; i<nelts ; i++ ) {
		if ( (i==0) || (i==nelts-1) ) soln[j%3][i] = 0.0	;
		else {
		    soln[j%3][i] = r*r*soln[(j-1)%3][i+1]	;
		    soln[j%3][i]+= 2.0*(1.0-r*r)*soln[(j-1)%3][i]	;
		    soln[j%3][i]+= r*r*soln[(j-1)%3][i-1]	;
		    soln[j%3][i]-= soln[(j-2)%3][i]	;
		}	// end of else clause
		System.out.printf("%8.4f ", soln[j%3][i]);
	    }	// end of for ( i ) loop
	    System.out.printf("\n")	;
	}	// end of for ( j ) loop
    }	// end of function main( )
}

//
//	EOF
//

