/*
 * Created on Feb 6, 2004
 *
 * To change the template for this generated file go to
 * Window - Preferences - Java - Code Generation - Code and Comments
 */
package pipe.modules.gspn;
//import java.applet.AudioClip;
//import java.applet.Applet.*;
import jama.Matrix;

import java.awt.Container;
import java.awt.event.ActionEvent;
import java.awt.event.ActionListener;
import java.text.DecimalFormat;
import java.util.ArrayList;

import javax.swing.BoxLayout;
import javax.swing.JDialog;

import pipe.dataLayer.Arc;
import pipe.dataLayer.DataLayer;
import pipe.dataLayer.PNMatrix;
import pipe.dataLayer.Place;
import pipe.dataLayer.PlaceTransitionObject;
import pipe.dataLayer.Transition;
import pipe.dataLayer.calculations.StateList;
import pipe.dataLayer.calculations.TreeTooBigException;
import pipe.dataLayer.calculations.myTree;
import pipe.gui.CreateGui;
import pipe.gui.widgets.ButtonBar;
import pipe.gui.widgets.PetriNetChooserPanel;
import pipe.gui.widgets.ResultsHTMLPane;
import pipe.modules.Module;
import pipe.modules.classification.Classification;


/**
 * @author Matthew Cook
 *
 * 
 */
public class GSPN extends Classification implements Module{
	
	private static final String MODULE_NAME = "GSPN Analysis";
	private PetriNetChooserPanel sourceFilePanel;
	private ResultsHTMLPane results;
  
//######################################################################################################################	
	public void run(DataLayer pnmlData) {
		// Build interface
		JDialog guiDialog = new JDialog(CreateGui.getApp(),MODULE_NAME,true);
		
		// 1 Set layout
		Container contentPane=guiDialog.getContentPane();
		contentPane.setLayout(new BoxLayout(contentPane,BoxLayout.PAGE_AXIS));
		
		// 2 Add file browser
		contentPane.add(sourceFilePanel=new PetriNetChooserPanel("Source net",pnmlData));
		
		// 3 Add results pane
		contentPane.add(results=new ResultsHTMLPane());
		
		// 4 Add button
		contentPane.add(new ButtonBar("Analyse GSPN",runAnalysis));
		//contentPane.add(new ButtonBar("Play a tune!",playMidi));
		// 5 Make window fit contents' preferred size
		guiDialog.pack();
		
		// 6 Move window to the middle of the screen
		guiDialog.setLocationRelativeTo(null);
		
		guiDialog.setVisible(true);
	}
//######################################################################################################################	
	/**
	 * Analyse button click handler
	 */
	ActionListener runAnalysis=new ActionListener() {
		public void actionPerformed(ActionEvent arg0) {
			DataLayer sourceDataLayer=sourceFilePanel.getDataLayer();
			String s="<h2>GSPN analysis results</h2>";
			if(sourceDataLayer==null) return;
			if (!hasTimedTransitions(sourceDataLayer)){
				s+= "This Petri net has no timed transitions, so GSPN analysis cannot be performed.";
				
			}
			else if (!hasImmediateTransitions(sourceDataLayer)){
				s+= "This Petri net has no immediate transitions, so GSPN analysis cannot be performed.";
				
			}
			
			
			else if (!testEqualConflict(sourceDataLayer)) {
				s+= "Condition Equal Conflict is not satisfied.  GSPN analysis cannot continue.";
				
			}
			
			else {
				
			//set up data for display	
			try {
				StateList reachSet = getReachabilitySet(sourceDataLayer);
				int reachSize = reachSet.size();
				getTransitionProbabilityMatrix(sourceDataLayer, reachSet);
				//int placeCount = sourceDataLayer.getPlaces().length;
				
				StateList vanishing = new StateList();
				
				StateList tangible = new StateList();
				getVanishingAndTangible (sourceDataLayer, reachSet, vanishing, tangible);
				if (tangible.size() == 0) {
					s+= "This petri net has no tangible states.  GSPN analysis cannot continue";
				}
				else  if (vanishing.size() == 0) {
					s+= "This petri net has no vanishing states.  GSPN analysis cannot continue";
				}
				else {
					double[][] c = new double [vanishing.size()][vanishing.size()];
					double[][] d = new double [vanishing.size()][tangible.size()];
					double[][] e = new double [tangible.size()][vanishing.size()];
					double[][] f = new double [tangible.size()][tangible.size()];
					c = probabilityMatrix(sourceDataLayer, vanishing, vanishing);
					d = probabilityMatrix(sourceDataLayer, vanishing, tangible);
					e = probabilityMatrix(sourceDataLayer, tangible, vanishing);
					f = probabilityMatrix(sourceDataLayer, tangible, tangible);
					Matrix cM = new Matrix(c);
					Matrix dM = new Matrix(d);
					Matrix eM = new Matrix(e);
					Matrix fM = new Matrix(f);
					int cSize = cM.getColumnDimension();
					Matrix iD = new Matrix (cSize, cSize);
					for (int i = 0; i<cSize; i++){
						for (int j = 0; j < cSize; j++){
							if (i==j) iD.set(i,j,1.0);
							else iD.set(i,j,0.0);
						}
					}
					Matrix iMinusCInverse = new Matrix(cSize, cSize);
					iMinusCInverse = iD.minus(cM);
					iMinusCInverse = iMinusCInverse.inverse();
					double[] sojournTimes = calcSojournTime(sourceDataLayer, tangible);
					double[] embeddedMarkovChainDist = getEmbeddedMarkovChainSteadyStateDistribution(cM,dM,eM,fM);
					Matrix meanVisits = calcMeanNumVisits(embeddedMarkovChainDist);
					double xHat = xHat(sojournTimes, embeddedMarkovChainDist);
					double[] meanCycleTimes = calcMeanCycleTimes(embeddedMarkovChainDist, xHat);
					int meanCycleLength = meanCycleTimes.length;
					double[] steadyStateDistribution = getSteadyStateDistribution(meanCycleTimes, sojournTimes);
					double[][] rates = rateMatrix(sourceDataLayer, tangible, vanishing);
					Matrix eTwiddle = new Matrix(rates);
					Matrix eTwiddleIMinusCInverse = eTwiddle.times(iMinusCInverse);
					Matrix steadyStateDistrib = new Matrix(1,meanCycleLength);
					for (int i=0; i< meanCycleLength; i++) {
						steadyStateDistrib.set(0,i,steadyStateDistribution[i]);
					}
					Matrix vanishingSteadyState = steadyStateDistrib.times(eTwiddleIMinusCInverse);
					double[] throughput = getTransitionThroughput(sourceDataLayer ,vanishing, tangible, vanishingSteadyState, steadyStateDistribution);
					
					s+=ResultsHTMLPane.makeTable(new String[]{
							"Entire reachability set",
							renderStateSpaceLinked(sourceDataLayer, reachSet)
						},1,false,false,true,false);
					
					s+=ResultsHTMLPane.makeTable(new String[]{
							"Set of vanishing markings",
							renderStateSpace(sourceDataLayer, vanishing)
					},1, false,false,true,false);
					
					s+=ResultsHTMLPane.makeTable(new String[]{
							"Set of tangible markings",
							renderStateSpace(sourceDataLayer, tangible)
					},1, false,false,true,false);
					
					s+=ResultsHTMLPane.makeTable(new String[]{
							"Probability of transition between two vanishing states",
							renderProbabilities(c, vanishing, vanishing)
					},1,false,false,true,false);	
										
					s+=ResultsHTMLPane.makeTable(new String[]{
							"Probability of transition from a vanishing state to a tangible state",
							renderProbabilities(d, vanishing, tangible)
					},1,false,false,true,false);
					
					s+=ResultsHTMLPane.makeTable(new String[]{
							"Probability of transition from a tangible state to a vanishing state",
							renderProbabilities(e, tangible, vanishing)
					},1,false,false,true,false);
					
					s+=ResultsHTMLPane.makeTable(new String[]{
							"Probability of transition from a tangible state to a tangible state",
							renderProbabilities(f, tangible, tangible)
					},1,false,false,true,false);
					
					s+=ResultsHTMLPane.makeTable(new String[]{
							"Steady state distribution  of the embedded Markov chain",
							renderLists(embeddedMarkovChainDist, tangible)
					},1,false,false,true,false);
					
					s+=ResultsHTMLPane.makeTable(new String[]{
							"Mean number of visits between tangible states",
							renderProbabilities(meanVisits.getArrayCopy(), tangible, tangible)
					},1,false,false,true,false);
					
					s+=ResultsHTMLPane.makeTable(new String[]{
							"Sojourn times for tangible states",
							renderLists(sojournTimes, tangible)
					},1,false,false,true,false);
					
					s+=ResultsHTMLPane.makeTable(new String[]{
							"Mean cycle times for tangible states",
							renderLists(meanCycleTimes, tangible)
					},1,false,false,true,false);
					
					s+=ResultsHTMLPane.makeTable(new String[]{
							"Steady state distribution for tangible states",
							renderLists(steadyStateDistribution, tangible)
					},1,false,false,true,false);			
					
					s+=ResultsHTMLPane.makeTable(new String[]{
							"Throughput",
							renderThroughput(sourceDataLayer, throughput)
					},1,false,false,true,false);			
				}
			} 	catch (TreeTooBigException e) {
				s+= e.getMessage();
				}
			}
			results.setText(s);
		}
		
	
		
			
		
	};
	/*
	ActionListener playMidi = new ActionListener() {

		public void actionPerformed(ActionEvent arg0) {
			// TODO Auto-generated method stub
			AudioClip midi = new AudioClip(CreateGui.imgPath+"FANTASY.mid");
		}
		
	};
	*/
//######################################################################################################################
	
	//Format StateList data nicely.
	private String renderStateSpace(DataLayer pnmldata, StateList data) {
		if((data.size()==0)||(data.get(0).length==0)) return "n/a";
		int markSize = data.get(0).length;
		ArrayList result=new ArrayList();
		// add headers to table
		result.add("");
		for (int i=0;i<markSize;i++) 
			{result.add(pnmldata.getPlaces()[i].getName());
			//result.add("<A NAME= 'M" + i + "'></A>");
			}
		
		for (int i=0; i<data.size(); i++) {
			result.add(data.getID(i));
			for (int j=0; j<markSize; j++)
				result.add(Integer.toString(data.get(i)[j]));
		}
		
		return ResultsHTMLPane.makeTable(result.toArray(),markSize + 1,false,true,true,true);
	}
	
	
	
//######################################################################################################################
	
	//Format StateList data nicely.
	private String renderStateSpaceLinked(DataLayer pnmldata, StateList data) {
		if((data.size()==0)||(data.get(0).length==0)) return "n/a";
		int markSize = data.get(0).length;
		ArrayList result=new ArrayList();
		// add headers to table
		result.add("");
		for (int i=0;i<markSize;i++) 
		{result.add(pnmldata.getPlaces()[i].getName());
		//result.add("<A NAME= 'M" + i + "'></A>");
		}
		
		for (int i=0; i<data.size(); i++) {
			result.add(data.getID(i)+ "<A NAME= 'M" + i + "'></A>");
			for (int j=0; j<markSize; j++)
				result.add(Integer.toString(data.get(i)[j]));
		}
		
		return ResultsHTMLPane.makeTable(result.toArray(),markSize + 1,false,true,true,true);
	}
	
//######################################################################################################################
	
	//Format throughput data nicely.
	private String renderThroughput(DataLayer pnmldata, double[] data) {
		if((data.length)==0) return "n/a";
		int transCount = data.length;
		ArrayList result=new ArrayList();
		// add headers to table
		result.add("Transition");
		result.add("Throughput");
    DecimalFormat f=new DecimalFormat();
    f.setMaximumFractionDigits(5);
		for (int i=0;i<transCount;i++) {
			result.add(pnmldata.getTransitions()[i].getName());
			result.add(f.format(data[i]));
		}
		
		return ResultsHTMLPane.makeTable(result.toArray(),2,false,true,true,true);
	}
	
//######################################################################################################################	
	
	//Format probability matrices nicely
	private String renderProbabilities(double[][] probabilities, StateList list1, StateList list2) {
		if((list1.size()==0)||(list2.get(0).length==0)) return "n/a";
		int rows = list1.size();
		int cols = list2.size();
		ArrayList result=new ArrayList();
		// add headers to table
		result.add("");
		for (int i=0;i<cols;i++) result.add("<A HREF='#M" + i + "'>" + list2.getID(i)+ "</A>");
		
    DecimalFormat f=new DecimalFormat();
    f.setMaximumFractionDigits(5);
		for (int i=0; i<rows; i++) {
			result.add("<A HREF='#M" + i + "'>" + list1.getID(i)+ "</A>");
			for (int j=0; j<cols; j++)
				result.add(f.format(probabilities[i][j]));
		}
		
		return ResultsHTMLPane.makeTable(result.toArray(),cols + 1,false,true,true,true);
	}
	
//######################################################################################################################	
	
	//Format lists of doubles nicely
	private String renderLists(double[] data, StateList list) {
		if((list.size()==0)) return "n/a";
		int rows = list.size();
		
		ArrayList result=new ArrayList();
		// add headers to table
		result.add("Marking");
		result.add("Value");

    DecimalFormat f=new DecimalFormat();
    f.setMaximumFractionDigits(5);
		for (int i=0; i<rows; i++) {
			result.add("<A HREF='#" + list.getID(i) + "'>" + list.getID(i).toString().toUpperCase()+ "</A>");
			result.add(f.format(data[i]));
		}
		
		return ResultsHTMLPane.makeTable(result.toArray(),2,false,true,true,true);
	}
	
//######################################################################################################################
	/**Qualitative analysis - see if the supplied GSPN is an EFC-GSPN.  This is a necessary precondition
	 * for quantitative analysis.
	 * @param DataLayer
	 * 
	 */
	private boolean isEFCGSPN (DataLayer pnmlData) {
		return extendedFreeChoiceNet(pnmlData);
	}
	
//######################################################################################################################	
	/**Generate the reachability set using myTree function
	 * Add each marking to an arraylist, testing to see if the
	 * marking is already present before adding.
	 * 
	 *
	 * @param DataLayer
	 * @return
	 */
	private StateList getReachabilitySet (DataLayer pnmlData) throws TreeTooBigException {
		int [][] fim = pnmlData.getForwardsIncidenceMatrix();
		int [][] bim = pnmlData.getBackwardsIncidenceMatrix();
		PNMatrix plus = new PNMatrix (fim);
		PNMatrix minus = new PNMatrix(bim);
		int[] marking = pnmlData.getCurrentMarkupMatrix();
		int markSize = pnmlData.getPlaces().length;
		StateList reachSetArray = new StateList();
		myTree reachSet = new myTree(marking, plus, minus, reachSetArray, pnmlData);
		return reachSetArray;
	}

//######################################################################################################################
	/**Get the initial marking of the supplied net
	 * 
	 * @param pnmlData
	 * @return
	 */
	private int[] getMarking(DataLayer pnmlData){
		
		int places = pnmlData.getPlaces().length;
		int[] marking = new int[places];
		for (int i = 0; i <places; i++){
			marking[i] = pnmlData.getPlace(i).getInitialMarking();
		}
		return marking;
	}
	
//######################################################################################################################	
	/**Caluculate whether a transition is enabled given a specific marking
	 * 
	 * @param DataLayer - the net  
	 * @param int[] - the marking
	 * @param int - the specific transition to test for enabled status
	 * @return boolean - an array of booleans specifying which transitions are enabled in the specified marking
	 */
	private boolean getTransitionEnabledStatus(DataLayer pnmlData, int[] marking, int transition) {
		int transCount = pnmlData.getTransitions().length;
		boolean[] result;
		result = new boolean[transCount];
		boolean answer;
		int[][] CMinus = pnmlData.getBackwardsIncidenceMatrix();
		int placeCount = pnmlData.getPlaces().length;
		
		for (int k = 0; k < transCount; k++) { //initialise matrix to true
			result[k] = true;
		}
		for (int i = 0; i <transCount;i++) {
			for (int j = 0; j <placeCount; j++){
				if (marking[j] < CMinus[j][i])
					result[i] = false;
			}
		}
		//print(result);
		return result[transition];
	}
	
//######################################################################################################################	
	/**Caluculate which transitions are enabled given a specific marking
	 * 
	 * @param DataLayer - the net  
	 * @param int[] - the marking
	 * @return boolean[] - an array of booleans specifying which transitions are enabled in the specified marking
	 */
	private boolean[] getTransitionEnabledStatusArray(DataLayer pnmlData, int[] marking) {
		int transCount = pnmlData.getTransitions().length;
		boolean[] result;
		result = new boolean[transCount];
		boolean hasTimed = false;
		boolean hasImmediate = false;
		int[][] CMinus = pnmlData.getBackwardsIncidenceMatrix();
		int placeCount = pnmlData.getPlaces().length;
		Transition[] transArray = pnmlData.getTransitions();
		
		for (int k = 0; k < transCount; k++) { //initialise matrix to true
			result[k] = true;
		}
		for (int i = 0; i <transCount;i++) {
			for (int j = 0; j <placeCount; j++){
				if (marking[j] < CMinus[j][i])
					result[i] = false;
			}
		}

		return result;
	}
//######################################################################################################################	
	private boolean[] getTangibleTransitionEnabledStatusArray(DataLayer pnmlData, int[] marking) {
		int transCount = pnmlData.getTransitions().length;
		boolean[] result;
		result = new boolean[transCount];
		boolean hasTimed = false;
		boolean hasImmediate = false;
		int[][] CMinus = pnmlData.getBackwardsIncidenceMatrix();
		int placeCount = pnmlData.getPlaces().length;
		Transition[] transArray = pnmlData.getTransitions();
		
		for (int k = 0; k < transCount; k++) { //initialise matrix to true
			result[k] = true;
		}
		for (int i = 0; i <transCount;i++) {
			for (int j = 0; j <placeCount; j++){
				if (marking[j] < CMinus[j][i])
					result[i] = false;
			}
		}
		for (int i = 0; i<transCount; i++) {
			if (transArray[i].getTimed()==true){
				hasTimed = true;
			}
			else {
				hasImmediate = true;
			}
		}
		if (hasTimed&&hasImmediate){
			for (int i = 0; i<transCount; i++){
				if (transArray[i].getTimed()==true){
					result[i] = false;
				}
			}
		}
		//print(result);
		return result;
	}
	
//######################################################################################################################
	/**Work out if a specified marking describes a tangible state.
	 * A state is either tangible (all enabled transitions are timed)
	 * or vanishing (there exists at least one enabled state that is transient, i.e. untimed).
	 * If an immediate transition exists, it will automatically fire before a timed transition.
	 * @param DataLayer - the net to be tested
	 * @param int[] - the marking of the net to be tested
	 * @return boolean - is it tangible or not
	 */

	private boolean isTangibleState(DataLayer pnmlData, int[] marking) {
		Transition[] trans = pnmlData.getTransitions();
		int numTrans = trans.length;
		boolean hasTimed = false;
		boolean hasImmediate = false;
		for (int i = 0; i < numTrans; i++ ){
			if (getTransitionEnabledStatus(pnmlData, marking, i) == true){
				if (trans[i].getTimed()== true){  
					//If any immediate transtions exist, the state is vanishing 
					//as they will fire immediately
				hasTimed = true;
				}
				else if (trans[i].getTimed()!= true) {
					hasImmediate = true;
				}
			}
		}
		if (hasTimed == true && hasImmediate == false)	
		return true;
		else
		return false;
	}
//######################################################################################################################
	/**This function takes a reachability set 
	 * and splits it into subsets of tangible and vanishing states
	 * @param DataLayer - the net to be processed
	 * @param StateList - the entire reachability set
	 * @param StateList - the list to be populated with vanishing states
	 * @param StateList - the list to be populated with tangible states
	 *
	 */
	private void getVanishingAndTangible(DataLayer pnmlData, StateList reachabilitySet, StateList vanishing, StateList tangible) {
		int size = reachabilitySet.size();
		for (int i = 0; i < size; i++) {
			String id = reachabilitySet.getID(i);
			if (isTangibleState(pnmlData, reachabilitySet.get(i))) 
				tangible.add(reachabilitySet.get(i),id);
			else
				vanishing.add(reachabilitySet.get(i),id);	
		}
	}
	

//######################################################################################################################
	/**This function calculates the mean number of visits at a particular transition, given a particular embedded Markov Process
	 * @param double[] - the embedded Markov Process
	 * @return Matrix
	 * 
	 */
	private Matrix calcMeanNumVisits (double[] embeddedMarkovProcess) {
		int size = embeddedMarkovProcess.length;
		Matrix result = new Matrix(size,size);
		for (int i = 0; i<size; i++) {
			for (int j = 0; j <size; j++){
			result.set(i,j,embeddedMarkovProcess[i]/embeddedMarkovProcess[j]);	
			}
		}
		return result;
	}
//######################################################################################################################
	/**This function determines the sojourn time for each state in a specified set of states.
	 * @param DataLater - the net to be analysed
	 * @param StateList - the list of tangible markings
	 * @return double[] - the array of sojourn times for each specific state
	 * 
	 */
	private double[] calcSojournTime (DataLayer pnmldata, StateList tangibleStates) {
		int numStates = tangibleStates.size();
		int numTrans = pnmldata.getTransitions().length;
		Transition[] trans = pnmldata.getTransitions();
		double[] sojournTime = new double[numStates];
	
		for (int i = 0; i < numStates; i++) {
		boolean[] transStatus = getTransitionEnabledStatusArray(pnmldata, tangibleStates.get(i));	
			double weights = 0;
			for (int j = 0; j <numTrans; j++) {
				if (transStatus[j] == true){
					weights += trans[j].getRate();
				}
			}
			sojournTime[i] = 1/weights;
		}
		return sojournTime;
	}
//######################################################################################################################
	
	
	//This is an intemediate calculation used to determine steady state distributions for tangible states.
	private double xHat(double[] piBar, double[] sojournTimes) {
		int size = piBar.length;
		double xHat = 0;
		for (int i = 0; i< size; i++) {
			xHat += (piBar[i] *sojournTimes[i]);
		}
		return xHat;
	}
	
//######################################################################################################################
	//This calculates the mean cycle times Tc(Mi) for tangible markings Mi
	private double[] calcMeanCycleTimes (double[] embeddedMarkovChain, double xHat) {
		int size = embeddedMarkovChain.length;
		double[] meanCycleTimes = new double[size];
		for (int i = 0; i < size; i++) {
			meanCycleTimes[i] = xHat/embeddedMarkovChain[i];
		}
		return meanCycleTimes;
		
	}
//######################################################################################################################
	private double[] getSteadyStateDistribution (double[] meanCycleTimes, double[] sojournTimes) {
		int size = meanCycleTimes.length;
		double[] steadyStateDistribution = new double[size];
		for (int i = 0; i<size; i++){
			steadyStateDistribution[i] = (sojournTimes[i]/meanCycleTimes[i]);
		}
		return steadyStateDistribution;
	}
	
//######################################################################################################################
	/**Test for condition Equal Conflict.  I.E., for all t1, t2
	 * in the set of transitions, where t1<>t2, that share the same
	 * input place, either t1, t2 are both in the set of timed transitions (T1)
	 * or t1, t2 are both in the set of immediate transitions (T2).
	 * 
	 * @param DataLayer
	 * @return boolean
	 */
	private boolean testEqualConflict (DataLayer pnmlData) {
		
		Place[] places = pnmlData.getPlaces();
		Arc[] arcs = pnmlData.getArcs();
		Transition[] trans = pnmlData.getTransitions();
		int arcsCount = arcs.length;
		int placesCount = places.length;
		
		for (int i = 0; i < placesCount ; i++)
		{
			boolean hasTimed = false;
			boolean hasUntimed = false;
			//get arcs with places[i] as source
			for (int j = 0; j < arcsCount; j++)
			{	
				if (arcs[j].getSource()==places[i]){
					PlaceTransitionObject targ = arcs[j].getTarget();
					if (((Transition)targ).getTimed() == true) {	
						hasTimed = true;
					} else {
					hasUntimed = true;
					}
				}
			if (hasTimed== true && hasUntimed == true)
				return false;
			}
		}
		return true;
	}
	
//######################################################################################################################
	/**Calculate the probability of changing from one marking to another
	 * Works out the intersection of transitions enabled to fire at a particular
	 * marking, transitions that can be reached from a particular marking and the 
	 * intersection of the two.  Then sums the firing rates of the intersection
	 * and divides it by the sum of the firing rates of the enabled transitions.
	 * @param DataLayer  
	 * @param int[] - marking1
	 * @param int[] - marking2
	 * @return double - the probability
	 */
	
	private double probMarkingAToMarkingB(DataLayer pnmlData, int[] marking1, int[] marking2){
	
		int markSize = marking1.length;
		int[][] incidenceMatrix = pnmlData.getIncidenceMatrix();
		int transCount = pnmlData.getTransitions().length;
		boolean[] marking1EnabledTransitions = new boolean[transCount];// = getTransitionEnabledStatus(pnmlData, marking1); //get list of transitions enabled at marking1
		boolean[] matchingTransition = new boolean[transCount];
		for (int i = 0; i <transCount; i++) {
		marking1EnabledTransitions[i] = getTransitionEnabledStatus(pnmlData, marking1, i);
		}
		
		//**************************************************** *************************************************
		for (int j = 0; j <transCount; j ++) {
			matchingTransition[j] = true;  //initialise matrix of potential transition values to true
		}
		//*****************************************************************************************************
		//get transition needed to fire to get from marking1 to marking2
		for (int i = 0; i < transCount; i++) {
			for (int k = 0; k <markSize; k++) {
			//if the sum of the incidence matrix and marking 1 doesn't equal marking 2, 
			//set that candidate transition possibility to be false 
				if (((int)marking1[k] + (int)incidenceMatrix[k][i])!= (int)marking2[k]){
					matchingTransition[i] = false;
					}
			}
		}
		//if the state is tangible, all transitions will be timed, 
		//so all can be considered in the probability calculation.
		//Otherwise, reset the enabled status of timed transitions to false, as immediate transitions
		//will always fire first.
		
		if (isTangibleState(pnmlData, marking1)== false) { 
			for (int i = 0; i <transCount; i++) {
				if (pnmlData.getTransitions()[i].getTimed() == true) {
					marking1EnabledTransitions[i] = false;
				}
			}
		}
		//*****************************************************************************************************
		//check if there are any potential transitions from marking 1 to marking 2 - if not return probability 0	
		boolean hasAChance = false;
		for (int i = 0; i < transCount; i++) {	
			if (matchingTransition[i] == true){
				hasAChance = true;	
			}
		}
		if (hasAChance == false) {
			return 0.0;
		}
		//*****************************************************************************************************
		//if a transition is needed to get to marking 2 but isn't enabled at marking 1, return probability 0
		boolean enabledAndMatching = false;
		for (int i = 0; i <transCount; i++) {	
			if ((matchingTransition[i] == true) && (marking1EnabledTransitions[i] == true)) {
				enabledAndMatching  = true;	
			}
		}
		if (enabledAndMatching = false) {
			return 0.0;
		}
		//******************************************************************************************************	
		//work out the sum of firing weights of input transitions
		double candidateTransitionWeighting = 0.0;
		for (int i = 0; i < transCount; i++) {
			if((matchingTransition[i] == true) && (marking1EnabledTransitions[i] == true)){
				candidateTransitionWeighting += pnmlData.getTransitions()[i].getRate(); 	
			}
		}
		//*****************************************************************************************************	
		//work out the sum of firing weights of enabled transitions
		double enabledTransitionWeighting = 0.0;
		for (int i = 0; i < transCount; i++) {
			if (marking1EnabledTransitions[i] == true) {
				enabledTransitionWeighting += pnmlData.getTransitions()[i].getRate();
			}
		}
		return (candidateTransitionWeighting/enabledTransitionWeighting);			
	}
		
//######################################################################################################################
	//This function is used to generate E-twiddle, the matrix specifying the rates for a specific state change leaving
	//a tangible Mi and entering a vanishing state Mj.
	private double getRateForSpecificStateChange(DataLayer pnmlData, int[] marking1, int[] marking2) {
		int markSize = marking1.length;
		int[][] incidenceMatrix = pnmlData.getIncidenceMatrix();
		int transCount = pnmlData.getTransitions().length;
		boolean[] marking1EnabledTransitions = new boolean[transCount];// = getTransitionEnabledStatus(pnmlData, marking1); //get list of transitions enabled at marking1
		boolean[] matchingTransition = new boolean[transCount];
		for (int i = 0; i <transCount; i++) {
			marking1EnabledTransitions[i] = getTransitionEnabledStatus(pnmlData, marking1, i);
		}
		
		//**************************************************** *************************************************
		for (int j = 0; j <transCount; j ++) {
			matchingTransition[j] = true;  //initialise matrix of potential transition values to true
		}
		//*****************************************************************************************************
		//get transition needed to fire to get from marking1 to marking2
		for (int i = 0; i < transCount; i++) {
			for (int k = 0; k <markSize; k++) {
				//if the sum of the incidence matrix and marking 1 doesn't equal marking 2, 
				//set that candidate transition possibility to be false 
				if (((int)marking1[k] + (int)incidenceMatrix[k][i])!= (int)marking2[k]){
					matchingTransition[i] = false;
				}
			}
		}
		//if the state is tangible, all transitions will be timed, 
		//so all can be considered in the probability calculation.
		//Otherwise, reset the enabled status of timed transitions to false, as immediate transitions
		//will always fire first.
		
		if (isTangibleState(pnmlData, marking1)== false) { 
			for (int i = 0; i <transCount; i++) {
				if (pnmlData.getTransitions()[i].getTimed() == true) {
					marking1EnabledTransitions[i] = false;
				}
			}
		}
		//*****************************************************************************************************
		//check if there are any potential transitions from marking 1 to marking 2 - if not return probability 0	
		boolean hasAChance = false;
		for (int i = 0; i < transCount; i++) {	
			if (matchingTransition[i] == true){
				hasAChance = true;	
			}
		}
		if (hasAChance == false) {
			return 0.0;
		}
		//*****************************************************************************************************
		//if a transition is needed to get to marking 2 but isn't enabled at marking 1, return probability 0
		boolean enabledAndMatching = false;
		for (int i = 0; i <transCount; i++) {	
			if ((matchingTransition[i] == true) && (marking1EnabledTransitions[i] == true)) {
				enabledAndMatching  = true;	
			}
		}
		if (enabledAndMatching = false) {
			return 0.0;
		}
		//******************************************************************************************************	
		//work out the sum of firing weights of input transitions
		double candidateTransitionWeighting = 0.0;
		for (int i = 0; i < transCount; i++) {
			if((matchingTransition[i] == true) && (marking1EnabledTransitions[i] == true)){
				candidateTransitionWeighting += pnmlData.getTransitions()[i].getRate(); 	
			}
		}
		
		return candidateTransitionWeighting;
		
		
	}
//######################################################################################################################
	//This function generates a matrix of e-twiddles - used in calculation of throughput.
	private double[][] rateMatrix(DataLayer pnmlData, StateList list1, StateList list2) {
		int rows = list1.size();
		int cols = list2.size();
		double[][] result = new double[rows][cols];
		for (int i = 0; i<rows; i++){
			for (int j = 0; j < cols; j++){
				result[i][j] = getRateForSpecificStateChange(pnmlData, list1.get(i),list2.get(j));
			}
		}
		return result;
	}
	
//######################################################################################################################
	//This function works out the throughput of an immediate transition for a vanishing state.
	private double getVanishingStateThroughput (DataLayer pnmldata, StateList list1, int transitionNumber, Matrix rateForSpecificState) {
		int length = list1.size();
		double result = 0;
		
		for (int i = 0; i< length; i++){
			double enabledTransitionRates = 0;
			double specifiedTransitionRate = 0;
			boolean[] transStatus = getTangibleTransitionEnabledStatusArray(pnmldata, list1.get(i));
			if (transStatus[transitionNumber]==true){
				int transCount = transStatus.length;
				for (int j = 0; j<transCount; j++){
					if (transStatus[j]==true){
						enabledTransitionRates += pnmldata.getTransitions()[j].getRate();
					}
				}
				specifiedTransitionRate = pnmldata.getTransitions()[transitionNumber].getRate();
				result+= (specifiedTransitionRate/enabledTransitionRates)*rateForSpecificState.get(0,i);
			}
			
		}
		return result;
	}
	
//######################################################################################################################
	
	private double getTransitionThroughputSPN (DataLayer pnmldata, StateList list, double[] steadyStateDistrib, int transitionNumber){
		int length = list.size();
		double result = 0;
		
		for (int i = 0; i< length; i++){
			double specifiedTransitionRate = 0;
			boolean[] transStatus = getTransitionEnabledStatusArray(pnmldata, list.get(i));
			//System.out.println(transStatus[0] + " " + transStatus[1]+ " "+ transStatus[2]+ " "+ transStatus[3]+ " "+ transStatus[4]+ " transStatuses" );
			if (transStatus[transitionNumber]==true){
				specifiedTransitionRate = pnmldata.getTransitions()[transitionNumber].getRate();
				//System.out.println(specifiedTransitionRate +" specified transition rate");
				result+= (specifiedTransitionRate*steadyStateDistrib[i]);
			}
		}
		return result;
	}
	
//######################################################################################################################	
	
	private double[] getTransitionThroughput(DataLayer pnmldata, StateList vanishing, StateList tangible, Matrix rateForSpecificChange, double[] steadyStateDistribution) {
		Transition[] transitions = pnmldata.getTransitions();
		int transCount = transitions.length;
		double[] result = new double[transCount];
		for (int i = 0; i<transCount; i++) {
			if (transitions[i].getTimed()==true){
				result[i] = getTransitionThroughputSPN(pnmldata, tangible, steadyStateDistribution, i);
			}
			else 
				result[i] = getVanishingStateThroughput(pnmldata, vanishing, i, rateForSpecificChange);
		}
		return result;
	}
	
	
//######################################################################################################################	
	/**Constructs a matrix of probabilities of changing from one marking to another.
	 * Uses the reachability set to determine all markings, and applies the 
	 * probMarkingAToMarkingB function to calculate probability for each pair of markings
	 * @param pnmlData
	 * @return
	 */
	private double[][]getTransitionProbabilityMatrix(DataLayer pnmlData, StateList reachabilitySet) {
		int setLength = reachabilitySet.size();
		double [][] transitionProbabilityMatrix = new double[setLength][setLength];
		int recordSize = reachabilitySet.get(0).length;
		int [] stateSpace1 = new int[recordSize];
		int [] stateSpace2 = new int [recordSize];
		for (int i = 0; i< setLength; i++){
			stateSpace1 = reachabilitySet.get(i);
			for (int j = 0; j < setLength; j++) {
				stateSpace2 = reachabilitySet.get(j);
				transitionProbabilityMatrix[i][j] = probMarkingAToMarkingB(pnmlData, stateSpace1, stateSpace2);
			}
		}
		return transitionProbabilityMatrix;
	}

//######################################################################################################################
	//This is a debugging function for viewing results in the console - not part of analysis 
	private void print (boolean[] transitions) {
		int size = transitions.length;
		for (int i = 0; i < size ; i++) {
			System.out.print( transitions[i] +" ");
		}
		System.out.println();
	}
	

//######################################################################################################################
	public String getName() {
		return MODULE_NAME;
	}

	
	
//######################################################################################################################
	/*
	private void runAnalysis(DataLayer pnmldata) {
		Transition[] transitions = pnmldata.getTransitions();
		int transCount = transitions.length;
		boolean hasTimed = false; */
		/*
		for (int i = 0; i< length; i++) {
			System.out.println("Timed attribute of transition " + transCount[i].getId());
			System.out.println(transCount[i].getTimed());
			if (transCount[i].getTimed()==true)
				hasTimed = true;
		}
		*/
/*		if (hasTimed = false){
			System.out.println("There are no timed transitions in this net.");
			System.out.println("GSPN analysis cannot be performed.");
			return;
		}
		
		if (testEqualConflict(pnmldata) == false) {
			System.out.println("Condition equal conflict does not hold, so analysis cannot continue");
			//TODO - put exit condition and apology in here
			return;
		}
		
		//if (isEFCGSPN(pnmldata) )
		System.out.println(isEFCGSPN(pnmldata) +" is EFCGSPN");
		System.out.println("Condition equal conflict = " + testEqualConflict(pnmldata));
		StateList reachSet = getReachabilitySet(pnmldata);
		int reachSize = reachSet.size();
		
		System.out.println("********************************************************");
		System.out.println("********************************************************");
		System.out.println("Decide whether states are vanishing or not");
		
		for (int j = 0; j <reachSize; j ++){
			int[] currentMarking = reachSet.get(j);
			System.out.println(isTangibleState(pnmldata, currentMarking) + "...is tangible state");
			printMarking(currentMarking);
		}
		
		getTransitionProbabilityMatrix(pnmldata, reachSet);
*/		
	
		//generate matrices C, D, E, F
		
/*		StateList tangible = new StateList();
		StateList vanishing = new StateList();
		
		getVanishingAndTangible(pnmldata, reachSet, vanishing, tangible);
*/		
		//MarkingProbability dummy = new MarkingProbability();
/*		double[][] cMatrix, dMatrix, eMatrix, fMatrix, iMinusC, pPrime;
		int cSize = vanishing.size();
		int fSize = tangible.size();
		
		cMatrix = probabilityMatrix(pnmldata, vanishing, vanishing);
		dMatrix = probabilityMatrix(pnmldata, vanishing, tangible);
		eMatrix = probabilityMatrix (pnmldata, tangible, vanishing);
		fMatrix = probabilityMatrix (pnmldata, tangible, tangible);
		iMinusC = identityMatrix(cSize);	//initialise to be the identity matrix
				

		System.out.println("C, D, E, F matrices follow:");
		printMatrix(cMatrix);
		printMatrix(dMatrix);
		printMatrix(eMatrix);
		printMatrix(fMatrix);
		printMatrix(iMinusC);		
		
		System.out.println("Generated C matrix" + cSize);
		
		//generate P' = F + ((E x (I - C) ^-1) x D)
		subtractMatrix(iMinusC, cMatrix);  //get the (I - C) term in the above expression
		
		Matrix iMin = new Matrix(iMinusC);
		Matrix f = new Matrix (fMatrix);
		Matrix d = new Matrix (dMatrix);
		Matrix e = new Matrix (eMatrix);
		Matrix c = new Matrix (cMatrix);
		
		
		Matrix iMinusCInverse = iMin.inverse();
		
		iMinusCInverse.print(8,5);
		
		Matrix iMinusCInverseD = new Matrix(iMinusCInverse.getRowDimension(),d.getColumnDimension());
		iMinusCInverseD = iMinusCInverse.times(d);
		
		iMinusCInverseD.print(8,5);
		
		Matrix eIMinusCInverseD = new Matrix(e.getRowDimension(),iMinusCInverseD.getColumnDimension());
		
		eIMinusCInverseD = e.times(iMinusCInverseD);
		
		eIMinusCInverseD.print(8,5);
		
		f.plusEquals(eIMinusCInverseD);
		
		System.out.println("The answer...");
		f.print(8,5);
		
		Matrix piBarM = new Matrix(1,f.getColumnDimension());
		for (int i = 0; i <f.getColumnDimension(); i++){
		piBarM.set(0,i,1);	
		}
		piBarM.print(8,5);
		
		Matrix ps = piBarM.times(f);
		
		System.out.println("See what this does");
		ps.print(8,5);
		
		Matrix inv = f.transpose();
		System.out.println("The inverse");
		inv.print(8,5);
		
		for (int j = 0; j <inv.getRowDimension(); j ++) {
			inv.set(j, j, (inv.get(j,j)) - 1);
		}
		inv.print (8,5);
		int row = inv.getColumnDimension();
		
		Matrix enh = new Matrix(row+1, row+1);
		
		for (int i = 0; i <= row; i++){
			enh.set(0, i, 1);
		}
		for (int i = 1; i<= row; i++) {
			for(int j = 0; j < row; j ++) {
				enh.set(i, j, inv.get(i-1,j));
			}
		}
		
		for (int i = 1; i <= row; i++) {
			enh.set(i, row, 0);
		}
		
		enh.print(8,5);
		*/
		
		//This generates the steady state distribution of the embedded Markov process.
	/*	double[] embeddedMarkovSteadyStateDistrib = reduction(enh);
		for (int i = 0; i < embeddedMarkovSteadyStateDistrib.length; i++) {
			System.out.println("Result " + i + " = " + embeddedMarkovSteadyStateDistrib[i]);
		}
		
		Matrix meanVisits = calcMeanNumVisits(embeddedMarkovSteadyStateDistrib);
		
		meanVisits.print(8,5);
		
		double[] sojournTime = calcSojournTime(pnmldata, tangible);
		
		int sLength = sojournTime.length;
		System.out.println("Sojourn times");
		for (int i = 0; i< sLength; i++) {
			System.out.print ("Sojourn time " + i + " " + sojournTime[i] + " ");
		}
		
		System.out.println();
		
		double xHat = xHat(sojournTime, embeddedMarkovSteadyStateDistrib);
		System.out.println(xHat +" xHat");

		double[] meanCycleTimes = calcMeanCycleTimes(embeddedMarkovSteadyStateDistrib, xHat);
		
		int meanCycleLength = meanCycleTimes.length;
		for (int i = 0; i < meanCycleLength; i++) {
			System.out.println(meanCycleTimes[i] + " meanCycleTimes["+ i +"]");
		}
		
		double[] steadyStateDistribution = getSteadyStateDistribution(meanCycleTimes, sojournTime);
		for (int i = 0; i<meanCycleLength; i++) {
			System.out.println(steadyStateDistribution[i] + " steady state distribution["+ i + "]");
		}
		
		double[][] rates = rateMatrix(pnmldata, tangible, vanishing);
		
		Matrix eTwiddle = new Matrix(rates);
		System.out.println("E twiddle");
		eTwiddle.print(5,8);
		
		Matrix eTwiddleIMinusCInverse = eTwiddle.times(iMinusCInverse);
		eTwiddleIMinusCInverse.print(5,8);
		Matrix steadyStateDistrib = new Matrix(1,meanCycleLength);
		for (int i=0; i< meanCycleLength; i++) {
			steadyStateDistrib.set(0,i,steadyStateDistribution[i]);
		}
		Matrix vanishingSteadyState = steadyStateDistrib.times(eTwiddleIMinusCInverse);
		vanishingSteadyState.print(5,8);
		*/
		
		/*
		double tp1 = getVanishingStateThroughput(pnmldata, vanishing, 0, vanishingSteadyState );
		System.out.println("Transition 1 throughput: " + tp1);
		double tp2 = getVanishingStateThroughput(pnmldata, vanishing, 1, vanishingSteadyState );
		System.out.println("Transition 2 throughput: " + tp2);
		double tp3 = getVanishingStateThroughput(pnmldata, vanishing, 2, vanishingSteadyState );
		
		System.out.println("Transition 3 throughput: " + tp3);
		
		double tp4 = getTransitionThroughputSPN(pnmldata, tangible, steadyStateDistribution, 4);
		System.out.println("Transition 4 throughput: " + tp4);
		*/
		//double[] throughput = getTransitionThroughput(pnmldata, vanishing, tangible, vanishingSteadyState, steadyStateDistribution);
		//for (int i = 0; i<transCount; i++){
		//	System.out.print( throughput[i] + " ");
		//}
		//next - calculate throughput for tangible states as per SPN methodology
		//calculate throughput for vanishing states .
		
		
		
	//}
	
//######################################################################################################################	
	
	/**
	 * 
	 * @param Matrix c where C is the matrix of transition probabilities from vanishing to vanishing states
	 * @param Matrix d where D is the matrix of transition probabilities from vanishing to tangible states
	 * @param Matrix e where E is the matrix of transition probabilities from tangible to vanishing states
	 * @param Matrix f where F is the matrix of transition probabilities from tangible to tangible states
	 * @return double[] - the steady state distributions for Markov chain embedded in tangible states
	 * This function produces a solution for the following equations:
	 * P' = F + E*((I-C)^-1)*D 
	 * &pi; bar * P' = &pi; bar
	 * &Sigma; (&pi; bar (i)) &forall; i == 1
	 * where &pi; bar is the vector of the embedded Markov chain steady state distribution.
	 * See Falko Bause - Stochastic Petri Nets - An Introduction to the Theory p181.
	 * 
	 */
	private double[] getEmbeddedMarkovChainSteadyStateDistribution(Matrix c, Matrix d, Matrix e, Matrix f){
		//Part one - generate P'  (= F + E*((I-C)^-1)*D)
		
		int cSize = c.getRowDimension();
		Matrix iMinusC  = new Matrix(cSize, cSize); //initialise as an identity matrix
		for (int i = 0; i<cSize; i++){
			for (int j = 0; j < cSize; j++) {
				if (i==j) iMinusC.set(i,j,1.0);
				else iMinusC.set(i,j,0);
			}
		}
		iMinusC.minusEquals(c);
		Matrix iMinusCInverse = iMinusC.inverse();
		Matrix iMinusCInverseD = new Matrix(iMinusCInverse.getRowDimension(),d.getColumnDimension());
		iMinusCInverseD = iMinusCInverse.times(d); 
		Matrix eIMinusCInverseD = new Matrix(e.getRowDimension(),iMinusCInverseD.getColumnDimension());
		eIMinusCInverseD = e.times(iMinusCInverseD);	
		f.plusEquals(eIMinusCInverseD);
		
		//Part two - rearrange P' and pi_bar to prepare them for Gaussian reduction
		Matrix piBarM = new Matrix(1,f.getColumnDimension());
		for (int i = 0; i <f.getColumnDimension(); i++){
			piBarM.set(0,i,1);	
		}
		
		Matrix inverseF = f.transpose();
		for (int j = 0; j <inverseF.getRowDimension(); j ++) {
			inverseF.set(j, j, (inverseF.get(j,j)) - 1);
		}
		int row = inverseF.getColumnDimension();
		
		Matrix solutionMatrix = new Matrix(row+1, row+1);
		
		for (int i = 0; i <= row; i++){
			solutionMatrix.set(0, i, 1);
		}
		for (int i = 1; i<= row; i++) {
			for(int j = 0; j < row; j ++) {
				solutionMatrix.set(i, j, inverseF.get(i-1,j));
			}
		}
		
		for (int i = 1; i <= row; i++) {
			solutionMatrix.set(i, row, 0);
		}
		double[] embeddedMarkovSteadyStateDistrib = reduction(solutionMatrix);
		return embeddedMarkovSteadyStateDistrib;
	
	}
	
//######################################################################################################################	
	
	/**
	 * See if the supplied net has any timed transitions.
	 * @param DataLayer
	 * @return boolean
	 * @author Matthew
	 *
	 */
	private boolean hasTimedTransitions(DataLayer pnmldata){
		Transition[] transitions = pnmldata.getTransitions();
		int transCount = transitions.length;
		boolean hasTimed = false;
		int length = transitions.length;
		
		 for (int i = 0; i< length; i++) {
		 if (transitions[i].getTimed()==true)
		 hasTimed = true;
		 }
		 
		if (hasTimed == true)
			return true;
		else 
			return false;
		
	}
	
//######################################################################################################################	
	
	/**
	 * See if the supplied net has any timed transitions.
	 * @param DataLayer
	 * @return boolean
	 * @author Matthew
	 *
	 */
	private boolean hasImmediateTransitions(DataLayer pnmldata){
		Transition[] transitions = pnmldata.getTransitions();
		int transCount = transitions.length;
		boolean hasImmediate = false;
		int length = transitions.length;
		
		for (int i = 0; i< length; i++) {
			if (transitions[i].getTimed()==false)
				hasImmediate = true;
		}
		
		if (hasImmediate == true)
			return true;
		else 
			return false;
		
	}
	
//######################################################################################################################
	/**This function performs a Gaussian reduction on a given Matrix, returning an array of values representing the solution.
	 * @param Matrix - the matrix of coefficients to be solved
	 * @return double[] - the array of solutions
	 * 
	 */
	private double[] reduction(Matrix input) {
		int row = input.getRowDimension();
		int col = input.getColumnDimension();
		double[] result = new double[col-1];
		//initialise results to 0 - have n-1 unknowns in n equations, so result can be 1 less than size of input matrix.
		for (int i = 0; i<col-1; i++){
			result[i] = 0;
		}
		//***********************************************************
		//First stage - reduce matrix of coefficients by substitution
		//***********************************************************
		boolean reducedThisRow = false;
		
		//Start - first row should have 1 as each coefficient.  Test if second row has 0 as coeffiecient - if so, move on and swap.
		for (int i = 0; i < row - 1; i++){
			for (int j = i + 1; j <row; j++){
				if ((input.get(j, i)== 0.0)&& reducedThisRow == false ){ //if the element is 0 and we haven't already  
					int k = j;											//reduced a row, search down the list till we find one, then swap it into the current position
					while ((input.get(k,i)== 0.0)&& k < row - 1 ) {
						k++;
						}
					if (k < row) {
						swapRows(input, j, k);
					} else {
						throw new ArithmeticException("Not enough parameters to calculate result");
					}
				} else if (input.get(j,i)!=0.0){					//reduce the row coeffecients by arithmetic substitution.
				double factor = ((input.get(i, i))/(input.get(j, i)));
				//System.out.println(factor + " Factor");
				//System.out.println(input.get(i,i) +" input.get(i,i)");
				//System.out.println(input.get(j,i) +" input.get(i,j)");
				multiplyRow(input,j,factor);
				subtractRow(input,i,j);
				//input.print(8,5);				//if the coefficient is 0 and we've already performed a reduction in this pass
				reducedThisRow = true;			//take no action and move onto the next
				}								
			
			}
			reducedThisRow = false;					
		}
	//************************************
	//next stage - backwards substitution.
	//************************************
		for (int i= row - 2; i>= 0; i--) {
			double backSub = 0;
			for (int j = i+1; j<row-1; j++) {	
					backSub = backSub + (result[j]*input.get(i,j));
				}
			result[i] = (input.get(i,row - 1) - backSub)/input.get(i,i);
			}
		return result;
	}
	
//######################################################################################################################	
	//Helper function for reduction function
	private void swapRows (Matrix input, int row1, int row2) {
		int col = input.getColumnDimension();
		double temp;
		for (int i = 0; i < col; i++) {
			temp = input.get(row1,i);
			input.set(row1,i,input.get(row2,i));
			input.set(row2,i,temp);
		}
	}
	
//######################################################################################################################	
	//Helper function for reduction function
	private void multiplyRow(Matrix input, int row, double factor) {
		int col = input.getColumnDimension();
		for (int i = 0; i <col ; i++){
			double newVal = (input.get(row,i))*factor;
			input.set(row, i,newVal);
		}
	}
	
//######################################################################################################################	
	//subtract the values of row1 from the values of row2 
	private void subtractRow(Matrix input, int row1, int row2) {
		int col = input.getColumnDimension();
		for (int i = 0; i <col; i++) {
			double r1 = input.get(row1,i);
			double r2 = input.get(row2,i);
			//System.out.println(r1 +" r1 " + r2 + " r2 " + i + " i ");
			input.set(row2,i,(r2 - r1));
		}
	}
	
//######################################################################################################################
	
	private double[][] probabilityMatrix(DataLayer pnmldata, StateList list1, StateList list2) {
		int list1Length = list1.size();
		int list2Length = list2.size();
		double[][] result = new double[list1Length][list2Length];
		
		for(int i = 0; i <list1Length; i++){
			for(int j = 0; j <list2Length; j++){	
				result[i][j] = probMarkingAToMarkingB(pnmldata, list1.get(i), list2.get(j));
				//System.out.println(result[i][j]+ "  probability " + i + " " + j);
			}
		}

		return result;
		
	}

//######################################################################################################################
	
	private void printMatrix(double[][] input) {
		int rows = input.length;
		int cols = input[0].length;
		System.out.println("Printing a matrix of "+ rows +" rows and " + cols +" columns.");
		for (int i = 0; i<rows; i++) {
			for (int j = 0; j < cols; j++) {
				System.out.print(input[i][j] + " ");
			}
			System.out.println("");
		}
	}
	
//######################################################################################################################
	
	private void printMatrix(int[][] input) {
		int rows = input.length;
		int cols = input[0].length;
		System.out.println("Printing a matrix of "+ rows +" rows and " + cols +" columns.");
		for (int i = 0; i<rows; i++) {
			for (int j = 0; j < cols; j++) {
				System.out.print(input[i][j] + " ");
			}
			System.out.println("");
		}
	}	

//######################################################################################################################
	
	private void printMarking (int[] marking) {
	int rows = marking.length;
	System.out.print("Marking as follows: ");
	for (int i = 0; i < rows; i++) {
	System.out.print(marking[i] + " ");
	}
	System.out.println();
	}


}