Blake Rollins’ code in Java

This is Blake’s code in Java: an amazing piece of work that will compute numerical solutions to initial value problems for first order differential equations \big\{ y'=f(x,y), y(x_0)=y_0 \big\} using both Euler and Improved Euler. The snippet of code that brings the displayed graphs on screen is also included.

Once compiled, the program will ask you first for the function f(x,y) on the left-hand side of the first order differential equation y'=f(x,y). The second input is the actual solution of the differential equation (or any other function you want to compare to the numerical solutions)—if you don’t want to submit any solution for comparison, simply type skip at that point. Both expressions should follow Maple notation, and so far only sine, cosines, tangents and exponentials are included in Blake’s libraries, but more functions can be easily added later.

The next inputs are, in this order, the value of x_0, the corresponding value of y_0, the time step h and finally the number of steps n.

The following exemplifies a session with Blake’s compiled program:

dy/dx = y
y = 3*e^x
x0 = 0
y0 = 3
h = 0.01
n = 5

Euler’s Approximation for equation: y
(0.00, 3)
(0.01, 3.030)
(0.02, 3.06030)
(0.03, 3.0909030)
(0.04, 3.121812030)
(0.05, 3.15303015030)

Improved Euler’s Approximation for equation: y
(0.00, 3)
(0.01, 3.03015)
(0.02, 3.0606030075)
(0.03, 3.091362067725375)
(0.04, 3.12243025650601501875)
(0.05, 3.1538106805839004696884375)

And the corresponding output graph looks like this:

The complete code can be read below. If you have any question or comment, please feel free to do so in the comments section below, or contact Blake directly at rollinbg at email dot sc dot edu.

Parser.java

Parser is a recursive function that grabs a supposedly correct mathematical expression (in Maple notation) and searches for operators within. Once all are found, then it solves each expression and combines them appropriately to obtain a well-formed function that java can interpret and use.

It is based on David Cole’s parsing tool.

import java.lang.Math;
import java.math.BigDecimal;;

public class Parser {
	private static final char[] validOperators = {'^','/','*','+','-'};
	private static final String[] validTrig = {"sin","cos","tan"};
	private static final BigDecimal e = BigDecimal.valueOf(2.71828182845904523536028747135266249775724709369995);

    private Parser(){
         
    }

    private static BigDecimal evaluate(String leftSide, char oper, String rightSide, BigDecimal x, BigDecimal y) throws IllegalArgumentException{
        BigDecimal total = BigDecimal.ZERO;
        BigDecimal leftResult = BigDecimal.ZERO;
        BigDecimal rightResult = BigDecimal.ZERO;

        //Solve left-side
        int operatorLoc = findOperatorLocation(leftSide);
        if( operatorLoc > 0 && operatorLoc < leftSide.length()-1 ){
        	leftResult = evaluate(leftSide.substring(0,operatorLoc), leftSide.charAt(operatorLoc), leftSide.substring(operatorLoc+1,leftSide.length()), x, y);
        }
        else{
        	try{
        		leftResult = BigDecimal.valueOf(Double.parseDouble(leftSide));
            }
            catch(Exception ex){
            	long multiplier = 1;
            	if(leftSide.charAt(0) == '-'){
            		multiplier = -1;
            		leftSide = leftSide.substring(1);
            	}
                if(leftSide.equalsIgnoreCase("x")){
                	leftResult = x;
                }
                else if(leftSide.equalsIgnoreCase("y")){
                	leftResult = y;
                }
                else if(leftSide.equalsIgnoreCase("e")){
                	leftResult = e;
                }
                else if(leftSide.indexOf('(') != 0){
                	String trig = leftSide.substring(0, leftSide.indexOf('('));
                	BigDecimal number = BigDecimal.valueOf(Double.parseDouble(leftSide.substring(leftSide.indexOf('(')+1, leftSide.indexOf(')'))));
                	if(trig.equalsIgnoreCase("sin")){
                		leftResult = BigDecimal.valueOf(Math.sin(number.doubleValue()));
                	}
                	else if(trig.equalsIgnoreCase("cos")){
                		leftResult = BigDecimal.valueOf(Math.cos(number.doubleValue()));
                	}
                	else if(trig.equalsIgnoreCase("tan")){
                		leftResult = BigDecimal.valueOf(Math.tan(number.doubleValue()));
                	}
                }
                else{
                	throw new IllegalArgumentException("Invalid value found in portion of equation: " + leftSide);
                }
                leftResult = leftResult.multiply(BigDecimal.valueOf(multiplier));
                //System.out.println("leftResult=" + leftResult);
            }
        }

        //Solve right-side
        operatorLoc = findOperatorLocation(rightSide);
        if( operatorLoc > 0 && operatorLoc < rightSide.length()-1 ){
            rightResult = evaluate(rightSide.substring(0,operatorLoc), rightSide.charAt(operatorLoc), rightSide.substring(operatorLoc+1,rightSide.length()), x, y);
        }
        else{
        	try{
        		rightResult = BigDecimal.valueOf(Double.parseDouble(rightSide));
            }
            catch(Exception ex){
            	long multiplier = 1;
            	if(rightSide.charAt(0) == '-'){
            		multiplier = -1;
            		rightSide = rightSide.substring(1);
            	}
            	if(rightSide.equalsIgnoreCase("x")){
                	rightResult = x;
                }
                else if(rightSide.equalsIgnoreCase("y")){
                	rightResult = y;
                }
                else if(rightSide.equalsIgnoreCase("e")){
                	rightResult = e;
                }
                else if(rightSide.indexOf('(') != 0){
                	String trig = rightSide.substring(0, rightSide.indexOf('('));
                	BigDecimal number = BigDecimal.valueOf(Double.parseDouble(rightSide.substring(rightSide.indexOf('(')+1, rightSide.indexOf(')'))));
                	if(trig.equalsIgnoreCase("sin")){
                		rightResult = BigDecimal.valueOf(Math.sin(number.doubleValue()));
                	}
                	else if(trig.equalsIgnoreCase("cos")){
                		rightResult = BigDecimal.valueOf(Math.cos(number.doubleValue()));
                	}
                	else if(trig.equalsIgnoreCase("tan")){
                		rightResult = BigDecimal.valueOf(Math.tan(number.doubleValue()));
                	}
                }
                else{
                	throw new IllegalArgumentException("Invalid value found in portion of equation: " + rightSide);
                }
            	rightResult = rightResult.multiply(BigDecimal.valueOf(multiplier));
            	//System.out.println("rightResult=" + rightResult);
            }
        }

        switch(oper)
        {
            case '/':
            	total = leftResult.divide(rightResult); break;
            case '*':
            	total = leftResult.multiply(rightResult); break;
            case '+':
            	total = leftResult.add(rightResult); break;
            case '-':
            	total = leftResult.subtract(rightResult); break;
            case '^':
            	total = BigDecimal.valueOf(Math.pow(leftResult.doubleValue(), rightResult.doubleValue())); break;
            default:
            	throw new IllegalArgumentException("Unknown operator.");
        }
        
        return total;
    }

    private static int findOperatorLocation(String string)
    {
        int index = -1;
        for(int i = validOperators.length-1; i >= 0; i--)
        {
            index = string.indexOf(validOperators[i]);
            if(index >= 0){
            	if(validOperators[i] == '-'){
            		boolean neg = false;
            		if(index == 0){
            			neg = true;
            		}
            		else{
            			for(int j = validOperators.length-1; j >= 0; j--){
                			if(string.charAt(index-1) == validOperators[j]){
                				neg = true;
                			}
                		}
            		}
            		if(!neg){
            			return index;
            		}
            		else{
            			string = string.substring(0, index) + "=" + string.substring(index+1);
            			return findOperatorLocation(string);
            		}
            	}
            	else{
            		return index;
            	}
            }
        }
        return index;
    }


    public static BigDecimal processEquation(String equation, BigDecimal x, BigDecimal y) throws IllegalArgumentException
    {
    	for(int i = 0; i < equation.length(); i++){
    		if(equation.charAt(i) == '('){
    			if(equation.indexOf(')') == 0){
    				throw new IllegalArgumentException("No close parentheses.");
    			}
    			else{
    				boolean flag = false;
    				if(equation.indexOf('(') > 2){
    					String trig = equation.substring(equation.indexOf('(')-3, equation.indexOf('('));
        				for(int j = 0; j < validTrig.length; j++){
        					if(trig.equalsIgnoreCase(validTrig[j])){
        						String left = equation.substring(0, equation.indexOf('(')+1);
        	    				String inner = equation.substring(equation.indexOf('(')+1, equation.indexOf(')'));
        	    				String right = equation.substring(equation.indexOf(')'));
        	    				equation = left + processEquation(inner, x, y) + right;
        	    				flag = true;
        					}
        				}
    				}
    				
    				if(!flag){
    					String left = equation.substring(0, equation.indexOf('('));
        				String inner = equation.substring(equation.indexOf('(')+1, equation.indexOf(')'));
        				String right = equation.substring(equation.indexOf(')')+1);
        				equation = left + processEquation(inner, x, y) + right;
    				}
    			}
    		}
    	}
    	return evaluate(equation,'+',"0", x, y);
    }
}

EulerMetod.java

import java.io.BufferedWriter;
import java.io.File;
import java.io.FileWriter;
import java.util.Scanner;
import javax.swing.*;
import java.util.ArrayList;
import java.math.BigDecimal;

public class EulerMethod {

	private static final int WIDTH = 500;
	private static final int HEIGHT = 500;
	private static BigDecimal h;
	private static BigDecimal x_0;
	private static BigDecimal y_0;
	private static BigDecimal n;
	private static String eq;
	private static String sol;
	private static ArrayList<ArrayList<Double>> coordinates;
	private static ArrayList<ArrayList<Double>> improvedCoordinates;
	private static ArrayList<ArrayList<Double>> solCoordinates;
	private static File dump;
	private static BufferedWriter out;
	
	public static void main(String[] args) {
		Scanner input = new Scanner(System.in);
		coordinates = new ArrayList<ArrayList<Double>>();
		improvedCoordinates = new ArrayList<ArrayList<Double>>();
		solCoordinates = new ArrayList<ArrayList<Double>>();
		System.out.println("Euler's Method:");
		System.out.print("dy/dx = ");
		eq = input.next();
		System.out.print("y = ");
		sol = input.next();
		if(sol.compareToIgnoreCase("skip") == 0){
			sol = "";
		}
		System.out.print("x0 = ");
		x_0 = input.nextBigDecimal();
		System.out.print("y0 = ");
		y_0 = input.nextBigDecimal();
		System.out.print("h = ");
		h = input.nextBigDecimal();
		System.out.print("n = ");
		n = input.nextBigDecimal();
		System.out.println("\nEuler's Approximation:");
		dump = new File("dump.txt");
		try{
			dump.createNewFile();
			FileWriter fw = new FileWriter("dump.txt");
			out = new BufferedWriter(fw);
			out.write("Euler's Approximation for equation: " + eq + "\n\n");
			getNextY(x_0.add((n.multiply(h))));
			out.write("\nImproved Euler's Approximation for equation: " + eq + "\n\n");
			System.out.println("\nImproved Euler's Approximation:");
			getNextYImproved(x_0.add((n.multiply(h))));
			if(sol.compareToIgnoreCase("") != 0){
				out.write("\nSolution for equation: " + sol + "\n\n");
				System.out.println("\nSolution:");
				solve(x_0);
			}
		}catch(Exception e){
			
		}
		Grapher grapher = new Grapher();
		grapher.setHeight(HEIGHT);
		grapher.setWidth(WIDTH);
		grapher.setCoordinates(coordinates);
		if(sol.compareTo("") != 0)
			grapher.setSolutionCoordinates(solCoordinates, sol);
		grapher.setImprovedCoordinates(improvedCoordinates);
		JFrame frame = new JFrame("Euler's Approximation");
		frame.setDefaultCloseOperation(JFrame.EXIT_ON_CLOSE);
		frame.add(grapher);
		frame.setSize(WIDTH+75, HEIGHT+75);
		frame.setVisible(true);
		try{
			out.close();
		}catch(Exception e){
			
		}
	}
	
	public static BigDecimal getNextY(BigDecimal x){
		ArrayList<Double> point = new ArrayList<Double>();
		
		if(x.compareTo(x_0) == 0){
			try{
				out.write("(" + x + ", " + y_0 + ")");
				out.newLine();
			}catch(Exception e){
				
			}
			System.out.println("(" + x + ", " + y_0 + ")");
			point.add(x.doubleValue());
			point.add(y_0.doubleValue());
			coordinates.add(point);
			return y_0;
		}
		else{
			BigDecimal y = getNextY(x.subtract(h));
			BigDecimal nextY = y.add((Parser.processEquation(eq, x.subtract(h), y).multiply(h)));
			//System.out.println("y=" + y + "\nnextY=" + nextY);
			if(nextY.doubleValue() == Double.NEGATIVE_INFINITY || nextY.doubleValue() == Double.POSITIVE_INFINITY){
				return nextY;
			}
			try{
				out.write("(" + x + ", " + nextY + ")");
				out.newLine();
			}catch(Exception e){
				
			}
			System.out.println("(" + x + ", " + nextY + ")");
			point.add(x.doubleValue());
			point.add(nextY.doubleValue());
			coordinates.add(point);
			return nextY;
		}
	}
	
	public static BigDecimal getNextYImproved(BigDecimal x){
		ArrayList<Double> point = new ArrayList<Double>();
		
		if(x.compareTo(x_0) == 0){
			try{
				out.write("(" + x + ", " + y_0 + ")");
				out.newLine();
			}catch(Exception e){
				
			}
			System.out.println("(" + x + ", " + y_0 + ")");
			point.add(x.doubleValue());
			point.add(y_0.doubleValue());
			improvedCoordinates.add(point);
			return y_0;
		}
		else{
			BigDecimal y = getNextYImproved(x.subtract(h));
			BigDecimal nextY = y.add((Parser.processEquation(eq, x.subtract(h), y).multiply(h)));
			BigDecimal nextYImproved = y.add((Parser.processEquation(eq, x.subtract(h), y).add(Parser.processEquation(eq, x, nextY)).divide(BigDecimal.valueOf(2)).multiply(h)));
			if(nextYImproved.doubleValue() == Double.NEGATIVE_INFINITY || nextYImproved.doubleValue() == Double.POSITIVE_INFINITY){
				return nextYImproved;
			}
			try{
				out.write("(" + x + ", " + nextYImproved + ")");
				out.newLine();
			}catch(Exception e){
				
			}
			System.out.println("(" + x + ", " + nextYImproved + ")");
			point.add(x.doubleValue());
			point.add(nextYImproved.doubleValue());
			improvedCoordinates.add(point);
			return nextYImproved;
		}
	}
	
	public static void solve(BigDecimal x){
		BigDecimal thisY;
		for(int i = 0; i <= n.intValue(); i++ ){
			thisY = Parser.processEquation(sol, x.add(h.multiply(BigDecimal.valueOf(i))), BigDecimal.ZERO);
			System.out.println("(" + x.add(h.multiply(BigDecimal.valueOf(i))) + ", " + thisY + ")");
			ArrayList<Double> point = new ArrayList<Double>();
			point.add(x.add(h.multiply(BigDecimal.valueOf(i))).doubleValue());
			point.add(thisY.doubleValue());
			solCoordinates.add(point);
		}
		
	}

}
  1. No comments yet.
  1. No trackbacks yet.

Leave a Reply

Fill in your details below or click an icon to log in:

WordPress.com Logo

You are commenting using your WordPress.com account. Log Out / Change )

Twitter picture

You are commenting using your Twitter account. Log Out / Change )

Facebook photo

You are commenting using your Facebook account. Log Out / Change )

Google+ photo

You are commenting using your Google+ account. Log Out / Change )

Connecting to %s

%d bloggers like this: