diff --git a/.DS_Store b/.DS_Store
new file mode 100644
index 00000000..221cb315
Binary files /dev/null and b/.DS_Store differ
diff --git a/.github/workflows/testrun.yml b/.github/workflows/testrun.yml
index 8adf1062..d40c41e6 100644
--- a/.github/workflows/testrun.yml
+++ b/.github/workflows/testrun.yml
@@ -20,13 +20,3 @@ jobs:
run: |
echo Running install
make install
- - name: "Run a simple graph"
- run: |
- cd scripts
- echo Trying to create graph
- ./creategraphfromdate.sh 2000 09
- - name: Upload graph image
- uses: actions/upload-artifact@v1
- with:
- name: Upload graph
- path: testrun/bview_2000_09/bview_2000_09_transparent.png
\ No newline at end of file
diff --git a/.gitignore b/.gitignore
index dee24d5d..6d917b17 100644
--- a/.gitignore
+++ b/.gitignore
@@ -1,5 +1,9 @@
tags
+tagsfile
*.o
+.*.un~
+.*.swp
+.swp
/bin/lglayout2D
/bin/lglayout3D
/bin/lglbreakup
@@ -24,4 +28,5 @@ src/lglbreakup
src/lglfileconvert
src/lglayout
src/lglrebuild
-**/*.pyx
\ No newline at end of file
+**/*.pyx
+.DS_Store
diff --git a/Java/Makefile b/Java/Makefile
index fad9a43c..0ef97a28 100644
--- a/Java/Makefile
+++ b/Java/Makefile
@@ -10,7 +10,7 @@ jar/%.jar: $(CLASS)
jarfiles: jar/ImageMaker.jar jar/LGLLib.jar jar/LGLView.jar
%.class : %.java
- javac -cp lib/Jama-1.0.3.jar:lib/VectorGraphics2D-0.13.jar src/Viewer2D/*.java src/ImageMaker/*.java
+ javac -cp lib/Jama-1.0.3.jar:lib/VectorGraphics2D-0.13.jar -cp src src/Viewer2D/*.java src/ImageMaker/*.java
all: $(CLASS)
diff --git a/Java/jar/ImageMaker.jar b/Java/jar/ImageMaker.jar
index 14efc7d8..e9fce429 100644
Binary files a/Java/jar/ImageMaker.jar and b/Java/jar/ImageMaker.jar differ
diff --git a/Java/jar/LGLLib.jar b/Java/jar/LGLLib.jar
index a59394bd..2951c849 100644
Binary files a/Java/jar/LGLLib.jar and b/Java/jar/LGLLib.jar differ
diff --git a/Java/jar/LGLView.jar b/Java/jar/LGLView.jar
index 77195ec0..c2d4c9bd 100644
Binary files a/Java/jar/LGLView.jar and b/Java/jar/LGLView.jar differ
diff --git a/Java/src/ImageMaker/GenerateImages.java b/Java/src/ImageMaker/GenerateImages.java
index 51315f57..56fc1eff 100755
--- a/Java/src/ImageMaker/GenerateImages.java
+++ b/Java/src/ImageMaker/GenerateImages.java
@@ -1,321 +1,335 @@
-//
-// Copyright (c) 2003 Alex Adai, All Rights Reserved.
-//
-// This program is free software; you can redistribute it and/or
-// modify it under the terms of the GNU General Public License as
-// published by the Free Software Foundation; either version 2 of
-// the License, or (at your option) any later version.
-//
-// This program is distributed in the hope that it will be useful,
-// but WITHOUT ANY WARRANTY; without even the implied warranty of
-// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
-// GNU General Public License for more details.
-//
-// You should have received a copy of the GNU General Public License
-// along with this program; if not, write to the Free Software
-// Foundation, Inc., 59 Temple Place, Suite 330, Boston,
-// MA 02111-1307 USA
-//
-
-package ImageMaker;
-
-import java.awt.Color;
-import java.awt.image.BufferedImage;
-import java.io.File;
-import java.io.FileOutputStream;
-import java.io.IOException;
-import java.text.MessageFormat;
-import java.util.ArrayList;
-import java.util.List;
-
-import Viewer2D.EdgesPanel;
-import Viewer2D.FormatVertex;
-import Viewer2D.ViewerIO;
-import de.erichseifert.vectorgraphics2d.Document;
-import de.erichseifert.vectorgraphics2d.Processor;
-import de.erichseifert.vectorgraphics2d.Processors;
-import de.erichseifert.vectorgraphics2d.VectorGraphics2D;
-import de.erichseifert.vectorgraphics2d.intermediate.CommandSequence;
-import de.erichseifert.vectorgraphics2d.util.PageSize;
-
-/**
- * SESS - 2014.05.10:
- *
- *
Cleanup code a little bit.
- *
Allow user to include edge color file ("-c <file>" at the end of
- * the command line).
- *
- *
- */
-public class GenerateImages {
- public static Color EDGE_COLOR = Color.white;
- public static String EDGE_COLOR_FILE = "color_file";
-
- public static void main(String[] args) {
- boolean loadedEdgeColors = false;
- if (args.length == 0) {
- message();
- }
-
- int[] windowSizes = new int[2];
- windowSizes[0] = new Integer(args[0]).intValue();
- windowSizes[1] = new Integer(args[1]).intValue();
- System.out.println("Loading flindeberg mod 2");
- System.out.println("Image size is " + windowSizes[0] + " x "
- + windowSizes[1]);
-
- // Check params
- String edgeFile = args[2];
- String edgeColorFile = EDGE_COLOR_FILE;
- List coordFiles = new ArrayList();
- boolean colorFileSwitch = false;
- for (int i = 3; i < args.length; i++) {
- String arg = args[i];
- if ("-c".equals(arg)) {
- colorFileSwitch = true;
- continue;
- }
- if (colorFileSwitch) {
- edgeColorFile = arg;
- break;
- }
- coordFiles.add(arg);
- }
- if (coordFiles.isEmpty()) {
- message();
- }
-
- System.out.println("Loading edge file: " + edgeFile + "...");
- ViewerIO verterIO = null;
- try {
- verterIO = new ViewerIO(new File(edgeFile));
- verterIO.loadSHORTFile();
- } catch (IOException e) {
- System.out.println(e.getMessage());
- System.exit(1);
- }
- System.out.println("Trying to load edge color file: " + edgeColorFile
- + "...");
- try {
- verterIO.loadEdgeColorFile(new File(edgeColorFile));
- loadedEdgeColors = true;
- } catch (java.io.FileNotFoundException e) {
- loadedEdgeColors = false;
- System.out.println(e.getMessage());
- } catch (java.io.IOException e) {
- System.out.println(e.getMessage());
- System.exit(1);
- }
- System.out.println("Edges loading complete.");
-
- generateDark(loadedEdgeColors, windowSizes, coordFiles, verterIO);
-
- generateLight(loadedEdgeColors, windowSizes, coordFiles, verterIO);
-
- generateTransparent(loadedEdgeColors, windowSizes, coordFiles, verterIO);
-
- // System.out.println("Going for vector graphics (currently broken)");
- // // Lets process coords files
- // for (String coordFile : coordFiles) {
- // try {
- // System.out.println("Loading " + coordFile + "...");
- // verterIO.loadVertexCoords(new File(coordFile));
-
-
- // FormatVertex formatter = new FormatVertex(
- // verterIO.getVertices(), verterIO.getStats(),
- // windowSizes, 1);
-
- // EdgesPanel panel = new EdgesPanel(verterIO.getEdges(),
- // verterIO.getVertices(), windowSizes[0], windowSizes[1]);
-
- // if (loadedEdgeColors)
- // panel.addEdgeColors(verterIO.getEdgeColorMap());
-
- // panel.showVertices(true);
- // panel.setVisibilityTest(true);
- // panel.setFormatter(formatter);
- // panel.setEdgeColor(EDGE_COLOR);
- // panel.setVertexColor(Color.white);
- // panel.setBackgroundColor(new Color(0f,0f,0f,0f));
-
- // // Use vector instead
- // VectorGraphics2D g2 = new VectorGraphics2D();
-
- // // Now the image has to be fitted to the given region
- // panel.fitData();
- // panel.writeVectorImage(g2);
-
- // CommandSequence commands = ((VectorGraphics2D) g2).getCommands();
-
- // String[] elements = {"pdf", "eps", "svg"};
- // for (String ext : elements) {
-
- // Processor pdfProcessor = Processors.get(ext);
- // Document doc = pdfProcessor.getDocument(commands, PageSize.A4);
-
- // String pngFile = MessageFormat.format(
- // "{0}_{1,number,0}x{2,number,0}_transparent.{3}", coordFile,
- // windowSizes[0], windowSizes[1], ext);
-
- // System.out.println("Preparing " + pngFile + "...");
- // try {
- // doc.writeTo(new FileOutputStream(pngFile));
- // } catch (Exception e) {
- // System.out.println("Could not write vector graphics");
- // }
-
- // System.out.println("Done.");
- // }
-
-
- // } catch (IOException e) {
- // System.out.println(MessageFormat.format(
- // "Error processing {0}:\n{1}", e.getMessage()));
- // }
- // }
- }
-
- private static void generateTransparent(boolean loadedEdgeColors, int[] windowSizes, List coordFiles,
- ViewerIO verterIO) {
- System.out.println("Going for transparent");
- // Lets process coords files
- for (String coordFile : coordFiles) {
- try {
- System.out.println("Loading " + coordFile + "...");
- verterIO.loadVertexCoords(new File(coordFile));
-
- String pngFile = MessageFormat.format(
- "{0}_{1,number,0}x{2,number,0}_transparent.png", coordFile,
- windowSizes[0], windowSizes[1]);
- System.out.println("Preparing " + pngFile + "...");
- FormatVertex formatter = new FormatVertex(
- verterIO.getVertices(), verterIO.getStats(),
- windowSizes, 1);
-
- EdgesPanel panel = new EdgesPanel(verterIO.getEdges(),
- verterIO.getVertices(), windowSizes[0], windowSizes[1]);
-
- if (loadedEdgeColors)
- panel.addEdgeColors(verterIO.getEdgeColorMap());
-
- panel.showVertices(true);
- panel.setVisibilityTest(true);
- panel.setFormatter(formatter);
- panel.setEdgeColor(EDGE_COLOR);
- panel.setVertexColor(Color.white);
- panel.setBackgroundColor(new Color(0f,0f,0f,0f));
-
- BufferedImage bufferedImage = new BufferedImage(windowSizes[0],
- windowSizes[1], BufferedImage.TYPE_INT_ARGB);
-
- // Now the image has to be fitted to the given region
- panel.fitData();
- panel.writeImage(pngFile, bufferedImage);
- System.out.println("Done.");
- } catch (IOException e) {
- System.out.println(MessageFormat.format(
- "Error processing {0}:\n{1}", e.getMessage()));
- }
- }
- }
-
- private static void generateLight(boolean loadedEdgeColors, int[] windowSizes, List coordFiles,
- ViewerIO verterIO) {
- System.out.println("Going for light");
- // Lets process coords files
- for (String coordFile : coordFiles) {
- try {
- System.out.println("Loading " + coordFile + "...");
- verterIO.loadVertexCoords(new File(coordFile));
-
- String pngFile = MessageFormat.format(
- "{0}_{1,number,0}x{2,number,0}_light.png", coordFile,
- windowSizes[0], windowSizes[1]);
- System.out.println("Preparing " + pngFile + "...");
- FormatVertex formatter = new FormatVertex(
- verterIO.getVertices(), verterIO.getStats(),
- windowSizes, 1);
-
- EdgesPanel panel = new EdgesPanel(verterIO.getEdges(),
- verterIO.getVertices(), windowSizes[0], windowSizes[1]);
-
- if (loadedEdgeColors)
- panel.addEdgeColors(verterIO.getEdgeColorMap());
-
- panel.showVertices(true);
- panel.setVisibilityTest(true);
- panel.setFormatter(formatter);
- panel.setEdgeColor(EDGE_COLOR);
- panel.setVertexColor(Color.white);
- panel.setBackgroundColor(Color.white);
-
- BufferedImage bufferedImage = new BufferedImage(windowSizes[0],
- windowSizes[1], BufferedImage.TYPE_INT_ARGB);
-
- // Now the image has to be fitted to the given region
- panel.fitData();
- panel.writeImage(pngFile, bufferedImage);
- System.out.println("Done.");
- } catch (IOException e) {
- System.out.println(MessageFormat.format(
- "Error processing {0}:\n{1}", e.getMessage()));
- }
- }
- }
-
- private static void generateDark(boolean loadedEdgeColors, int[] windowSizes, List coordFiles,
- ViewerIO verterIO) {
- System.out.println("Going for dark");
- // Lets process coords files
- for (String coordFile : coordFiles) {
- try {
- System.out.println("Loading " + coordFile + "...");
- verterIO.loadVertexCoords(new File(coordFile));
-
- String pngFile = MessageFormat.format(
- "{0}_{1,number,0}x{2,number,0}_dark.png", coordFile,
- windowSizes[0], windowSizes[1]);
- System.out.println("Preparing " + pngFile + "...");
- FormatVertex formatter = new FormatVertex(
- verterIO.getVertices(), verterIO.getStats(),
- windowSizes, 1);
-
- EdgesPanel panel = new EdgesPanel(verterIO.getEdges(),
- verterIO.getVertices(), windowSizes[0], windowSizes[1]);
-
- if (loadedEdgeColors)
- panel.addEdgeColors(verterIO.getEdgeColorMap());
-
- panel.showVertices(true);
- panel.setVisibilityTest(true);
- panel.setFormatter(formatter);
- panel.setEdgeColor(EDGE_COLOR);
- panel.setVertexColor(Color.white);
- panel.setBackgroundColor(Color.BLACK);
-
- BufferedImage bufferedImage = new BufferedImage(windowSizes[0],
- windowSizes[1], BufferedImage.TYPE_INT_ARGB);
-
- // Now the image has to be fitted to the given region
- panel.fitData();
- panel.writeImage(pngFile, bufferedImage);
- System.out.println("Done.");
- } catch (IOException e) {
- System.out.println(MessageFormat.format(
- "Error processing {0}:\n{1}", e.getMessage()));
- }
- }
- }
-
- public static void message() {
- System.out
- .println("Arguments:\n\n"
- + "\t... [-c ]\n\n"
- + "If no colors file specified program will try to load file named \""
- + EDGE_COLOR_FILE + "\".\n"
- + "By default edges are white. flindeberg mod");
- System.exit(1);
- }
-
-}
+package ImageMaker;
+
+import java.awt.Color;
+import java.awt.image.BufferedImage;
+import java.io.File;
+import java.io.FileNotFoundException;
+import java.io.IOException;
+import java.text.MessageFormat;
+import java.util.ArrayList;
+import java.util.List;
+
+import Viewer2D.EdgesPanel;
+import Viewer2D.FormatVertex;
+import Viewer2D.ViewerIO;
+
+public class GenerateImages {
+
+ public static boolean GENERATE_IN_SEPARATE_THREAD = true;
+
+ /**
+ * Entry point.
+ *
+ * @param args
+ */
+ public static void main(String[] args) {
+ ParseArguments pa = new ParseArguments(false);
+ pa.parse(args);
+
+ printStats(pa);
+
+ GENERATE_IN_SEPARATE_THREAD = true;
+
+ generate(pa, "dark with labels", "dark_withlabels", Color.BLACK,
+ true /*scale*/, true /*printLabels*/, true /*useAlignmentCenterArg*/);
+ generate(pa, "dark without labels", "dark_nolabels", Color.BLACK,
+ true /*scale*/, false /*printLabels*/, true /*useAlignmentCenterArg*/);
+
+ generate(pa, "dark without scale without labels", "dark_withoutscale_withoutlabels", Color.BLACK,
+ false /*scale*/, false /*printLabels*/, false /*useAlignmentCenterArg*/);
+ generate(pa, "dark without scale with labels", "dark_withoutscale_withlabels", Color.BLACK,
+ false /*scale*/, true /*printLabels*/, false /*useAlignmentCenterArg*/);
+
+ generate(pa, "light without scale without labels", "light_withoutscale_withoutlabels", Color.WHITE,
+ false /*scale*/, false /*printLabels*/, false /*useAlignmentCenterArg*/);
+ generate(pa, "transparent without scale without labels", "transparent_withoutscale_withoutlabels", new Color(0f, 0f, 0f, 0f),
+ false /*scale*/, false /*printLabels*/, false /*useAlignmentCenterArg*/);
+ }
+
+ private static void printStats(ParseArguments pa) {
+ try {
+ ViewerIO v = createViewerIO(pa);
+ v.loadVertexCoords(new File(pa.coordFiles.get(0)));
+ v.getStats().print();
+ } catch (Exception e) {
+ e.printStackTrace();
+ }
+ }
+
+ private static ViewerIO createViewerIO(ParseArguments pa) {
+ ViewerIO verterIO = null;
+ try {
+ verterIO = new ViewerIO(new File(pa.edgeFile));
+ verterIO.loadSHORTFile();
+ } catch (IOException e) {
+ System.out.println(e.getMessage());
+ System.exit(1);
+ }
+ verterIO.setLabelScale(pa.scaling);
+
+ verterIO.setMinMaxXY(pa.minX, pa.minY, pa.maxX, pa.maxY);
+
+ try {
+ verterIO.loadEdgeColorFile(new File(pa.edgeColorFile));
+ } catch (java.io.FileNotFoundException e) {
+ System.out.println(e.getMessage());
+ } catch (java.io.IOException e) {
+ System.out.println(e.getMessage());
+ System.exit(1);
+ }
+
+ return verterIO;
+ }
+
+ private static void generate(ParseArguments pa, String displayname, String name, Color background,
+ boolean scale, boolean printLabels, boolean useAlignmentCenterArg) {
+
+ ViewerIO verterIO = createViewerIO(pa);
+
+ Runnable task = new Runnable() {
+ @Override
+ public void run() {
+ generate(displayname, name, background, !verterIO.getEdgeColorMap().isEmpty(), pa.windowSizes, pa.coordFiles,
+ printLabels ? pa.labelFile : "", useAlignmentCenterArg ? pa.alignmentCenter : false, verterIO, scale);
+ }
+ };
+
+ if (GENERATE_IN_SEPARATE_THREAD) {
+ new Thread(task).start();
+ }
+ else {
+ task.run();
+ }
+ }
+
+ private static void generate(String displayname, String name, Color background, boolean loadedEdgeColors,
+ int[] windowSizes, List coordFiles, String labelFile, boolean alignmentCenter, ViewerIO verterIO, boolean scale) {
+ if (!labelFile.isEmpty())
+ loadLabels(labelFile, verterIO);
+ else
+ verterIO.clearLabels();
+
+ for (String coordFile : coordFiles) {
+ try {
+ verterIO.loadVertexCoords(new File(coordFile));
+
+ FormatVertex formatter = new FormatVertex(verterIO.getVertices(), verterIO.getLabels(),
+ verterIO.getLabelScale(), scale ? verterIO.getMinX() : 0, scale ? verterIO.getMinY() : 0,
+ scale ? verterIO.getMaxX() : 0, scale ? verterIO.getMaxY() : 0, alignmentCenter,
+ verterIO.getStats(), windowSizes, 1);
+
+ EdgesPanel panel = new EdgesPanel(verterIO.getEdges(), verterIO.getVertices(), verterIO.getLabels(),
+ windowSizes[0], windowSizes[1]);
+
+ if (loadedEdgeColors)
+ panel.addEdgeColors(verterIO.getEdgeColorMap());
+
+ panel.showVertices(true);
+ panel.setVisibilityTest(true);
+ panel.setFormatter(formatter);
+ panel.setEdgeColor(EDGE_COLOR);
+ panel.setVertexColor(Color.white);
+ panel.setBackgroundColor(background);
+
+ BufferedImage bufferedImage = new BufferedImage(windowSizes[0], windowSizes[1],
+ BufferedImage.TYPE_INT_ARGB);
+
+ // Now the image has to be fitted to the given region
+ panel.fitData();
+ String pngFile = MessageFormat.format("{0}_{1,number,0}x{2,number,0}_" + name + ".png", coordFile,
+ windowSizes[0], windowSizes[1]);
+ panel.writeImage(pngFile, bufferedImage);
+ System.out.println("Done.");
+ } catch (IOException e) {
+ System.out.println(MessageFormat.format("Error processing {0}:\n{1}", e.getMessage()));
+ }
+ }
+ }
+
+ static void loadLabels(String labelFile, ViewerIO verterIO) {
+ if (labelFile.isEmpty()) return;
+ System.out.println("Loading label file: " + labelFile + "...");
+ try {
+ verterIO.loadLabelFile(new File(labelFile));
+ } catch (FileNotFoundException e) {
+ System.out.println(e.getMessage());
+ } catch (IOException e) {
+ System.out.println(e.getMessage());
+ System.exit(1);
+ }
+ System.out.println("Labels loading complete.");
+ }
+
+ public static class ParseArguments {
+
+ public String edgeFile;
+ public String edgeColorFile;
+ public String labelFile;
+ public double scaling;
+ public double minX;
+ public double maxX;
+ public double minY;
+ public double maxY;
+ int[] windowSizes;
+ boolean alignmentCenter;
+ public List coordFiles;
+ private boolean viewer2d;
+
+ public ParseArguments(boolean viewer2d) {
+ this.viewer2d = viewer2d;
+
+ }
+
+ void usage() {
+ if (viewer2d)
+ message2();
+ else
+ message();
+
+ }
+
+ public void parse(String[] args) {
+ int minArguments = 3;
+
+ if (viewer2d) {
+ minArguments = 1;
+ if (args.length == 0)
+ return;
+ }
+
+ if (args.length < minArguments) {
+ usage();
+ }
+ int argno = 0;
+ if (!viewer2d) {
+ windowSizes = new int[2];
+ windowSizes[0] = Integer.parseInt(args[0]);
+ windowSizes[1] = Integer.parseInt(args[1]);
+ argno += 2;
+ System.out.println("Image size is " + windowSizes[0] + " x " + windowSizes[1]);
+ }
+
+ // Check params
+ edgeFile = args[argno++];
+ edgeColorFile = "";
+ coordFiles = new ArrayList();
+ labelFile = "";
+ boolean colorFileSwitch = false;
+ boolean labelFileSwitch = false;
+ boolean scaleSwitch = false;
+ scaling = 1;
+
+ boolean wasMin = false;
+ boolean wasMax = false;
+ boolean minSwitch = false;
+ boolean maxSwitch = false;
+ minX = 0;
+ maxX = 0;
+ minY = 0;
+ maxY = 0;
+
+ alignmentCenter = false;
+ boolean alignSwitch = false;
+ for (int i = argno; i < args.length; i++) {
+ String arg = args[i];
+ if ("-c".equals(arg)) {
+ colorFileSwitch = true;
+ continue;
+ }
+ if ("-l".equals(arg)) {
+ labelFileSwitch = true;
+ continue;
+ }
+ if ("-s".equals(arg)) {
+ scaleSwitch = true;
+ continue;
+ }
+ if ("-m".equals(arg) && !viewer2d) {
+ minSwitch = true;
+ continue;
+ }
+ if ("-a".equals(arg) && !viewer2d) {
+ alignSwitch = true;
+ continue;
+ }
+ if ("-M".equals(arg) && !viewer2d) {
+ maxSwitch = true;
+ continue;
+ }
+ if (scaleSwitch) {
+ scaleSwitch = false;
+ scaling = Double.parseDouble(arg);
+ continue;
+ }
+ if (minSwitch) {
+ minSwitch = false;
+ wasMin = true;
+ String[] a = arg.split(",");
+ if (a.length != 2) {
+ System.out.println("Error:-m requires exactly 2 coordinates");
+ System.exit(1);
+
+ }
+ minX = Double.parseDouble(a[0]);
+ minY = Double.parseDouble(a[1]);
+ continue;
+ }
+ if (maxSwitch) {
+ maxSwitch = false;
+ wasMax = true;
+ String[] a = arg.split(",");
+ if (a.length != 2) {
+ System.out.println("Error:-M requires exactly 2 coordinates");
+ System.exit(1);
+ }
+ maxX = Double.parseDouble(a[0]);
+ maxY = Double.parseDouble(a[1]);
+ continue;
+ }
+ if (alignSwitch) {
+ alignSwitch = false;
+ if (arg.equals("center"))
+ alignmentCenter = true;
+ continue;
+ }
+ if (labelFileSwitch) {
+ labelFile = arg;
+ labelFileSwitch = false;
+
+ continue;
+ }
+ if (colorFileSwitch) {
+ edgeColorFile = arg;
+ colorFileSwitch = false;
+ continue;
+ }
+ coordFiles.add(arg);
+ }
+ if (coordFiles.isEmpty()) {
+ usage();
+ }
+
+ if (wasMax ^ wasMin) {
+ System.out.println("Error:Both -m and -M need to be used at the same time, one of them is missing");
+ System.exit(1);
+ }
+ }
+ }
+
+ public static Color EDGE_COLOR = Color.white;
+ public static String EDGE_COLOR_FILE = "color_file";
+
+ public static void message() {
+ System.out.println("Arguments:\n\n"
+ + "\t... [-c ] [-l ] [-m minx,miny -M -maxx,maxy] [-a center]\n\n"
+ + "If no colors file specified program will try to load file named \"" + EDGE_COLOR_FILE + "\".\n"
+ + "By default edges are white. flindeberg mod");
+ System.exit(1);
+ }
+
+ public static void message2() {
+ System.out.println("Arguments:\n\n" + "\t [-c ] [-l ]\n\n"
+ + "If no colors file specified program will try to load file named \"" + EDGE_COLOR_FILE + "\".\n");
+ System.exit(1);
+ }
+}
\ No newline at end of file
diff --git a/Java/src/Jama/ChangeLog b/Java/src/Jama/ChangeLog
new file mode 100644
index 00000000..1557c282
--- /dev/null
+++ b/Java/src/Jama/ChangeLog
@@ -0,0 +1,42 @@
+2012-11-09 Bruce Miller
+ * Version 1.0.3 Released.
+ * Patched hqr2 method in Jama.EigenvalueDecomposition
+ to avoid infinite loop;
+ Thanks Frederic Devernay
+ * Updated unchecked type coding in Matrix read.
+ * Added serialization ID to Jama.Matrix and Decomposition classes.
+ * Cleaned up javadoc documentation to Jama.Matrix and Decomposition classes.
+2005-07-13 Bruce Miller
+ * Version 1.0.2 Released.
+
+2005-07-12 Bruce Miller
+ * Although SVD sometimes fails on cases m < n,
+ the test code, apparently successfully, invokes
+ SVD to compute rank on an m
+ * Version 1.0.1 Released.
+
+2000-09-11 Bruce Miller
+ * Jama.Matrix print methods which create a NumberFormat, now set its
+ Locale to US so that the reader will recognize them even when the
+ default locale is not US. Similar change to Jama.test.TestMatrix.
+ (Thanks Ulrich Eberhardinger )
+
+1998-08-05 The Jama Team
+ * Initial Version released (1.0.0)
+
diff --git a/Java/src/Jama/CholeskyDecomposition.java b/Java/src/Jama/CholeskyDecomposition.java
new file mode 100644
index 00000000..57f4ccfb
--- /dev/null
+++ b/Java/src/Jama/CholeskyDecomposition.java
@@ -0,0 +1,201 @@
+package Jama;
+
+ /** Cholesky Decomposition.
+
+ For a symmetric, positive definite matrix A, the Cholesky decomposition
+ is an lower triangular matrix L so that A = L*L'.
+
+ If the matrix is not symmetric or positive definite, the constructor
+ returns a partial decomposition and sets an internal flag that may
+ be queried by the isSPD() method.
+ */
+
+public class CholeskyDecomposition implements java.io.Serializable {
+
+/* ------------------------
+ Class variables
+ * ------------------------ */
+
+ /** Array for internal storage of decomposition.
+ @serial internal array storage.
+ */
+ private double[][] L;
+
+ /** Row and column dimension (square matrix).
+ @serial matrix dimension.
+ */
+ private int n;
+
+ /** Symmetric and positive definite flag.
+ @serial is symmetric and positive definite flag.
+ */
+ private boolean isspd;
+
+/* ------------------------
+ Constructor
+ * ------------------------ */
+
+ /** Cholesky algorithm for symmetric and positive definite matrix.
+ Structure to access L and isspd flag.
+ @param Arg Square, symmetric matrix.
+ */
+
+ public CholeskyDecomposition (Matrix Arg) {
+
+
+ // Initialize.
+ double[][] A = Arg.getArray();
+ n = Arg.getRowDimension();
+ L = new double[n][n];
+ isspd = (Arg.getColumnDimension() == n);
+ // Main loop.
+ for (int j = 0; j < n; j++) {
+ double[] Lrowj = L[j];
+ double d = 0.0;
+ for (int k = 0; k < j; k++) {
+ double[] Lrowk = L[k];
+ double s = 0.0;
+ for (int i = 0; i < k; i++) {
+ s += Lrowk[i]*Lrowj[i];
+ }
+ Lrowj[k] = s = (A[j][k] - s)/L[k][k];
+ d = d + s*s;
+ isspd = isspd & (A[k][j] == A[j][k]);
+ }
+ d = A[j][j] - d;
+ isspd = isspd & (d > 0.0);
+ L[j][j] = Math.sqrt(Math.max(d,0.0));
+ for (int k = j+1; k < n; k++) {
+ L[j][k] = 0.0;
+ }
+ }
+ }
+
+/* ------------------------
+ Temporary, experimental code.
+ * ------------------------ *\
+
+ \** Right Triangular Cholesky Decomposition.
+
+ For a symmetric, positive definite matrix A, the Right Cholesky
+ decomposition is an upper triangular matrix R so that A = R'*R.
+ This constructor computes R with the Fortran inspired column oriented
+ algorithm used in LINPACK and MATLAB. In Java, we suspect a row oriented,
+ lower triangular decomposition is faster. We have temporarily included
+ this constructor here until timing experiments confirm this suspicion.
+ *\
+
+ \** Array for internal storage of right triangular decomposition. **\
+ private transient double[][] R;
+
+ \** Cholesky algorithm for symmetric and positive definite matrix.
+ @param A Square, symmetric matrix.
+ @param rightflag Actual value ignored.
+ @return Structure to access R and isspd flag.
+ *\
+
+ public CholeskyDecomposition (Matrix Arg, int rightflag) {
+ // Initialize.
+ double[][] A = Arg.getArray();
+ n = Arg.getColumnDimension();
+ R = new double[n][n];
+ isspd = (Arg.getColumnDimension() == n);
+ // Main loop.
+ for (int j = 0; j < n; j++) {
+ double d = 0.0;
+ for (int k = 0; k < j; k++) {
+ double s = A[k][j];
+ for (int i = 0; i < k; i++) {
+ s = s - R[i][k]*R[i][j];
+ }
+ R[k][j] = s = s/R[k][k];
+ d = d + s*s;
+ isspd = isspd & (A[k][j] == A[j][k]);
+ }
+ d = A[j][j] - d;
+ isspd = isspd & (d > 0.0);
+ R[j][j] = Math.sqrt(Math.max(d,0.0));
+ for (int k = j+1; k < n; k++) {
+ R[k][j] = 0.0;
+ }
+ }
+ }
+
+ \** Return upper triangular factor.
+ @return R
+ *\
+
+ public Matrix getR () {
+ return new Matrix(R,n,n);
+ }
+
+\* ------------------------
+ End of temporary code.
+ * ------------------------ */
+
+/* ------------------------
+ Public Methods
+ * ------------------------ */
+
+ /** Is the matrix symmetric and positive definite?
+ @return true if A is symmetric and positive definite.
+ */
+
+ public boolean isSPD () {
+ return isspd;
+ }
+
+ /** Return triangular factor.
+ @return L
+ */
+
+ public Matrix getL () {
+ return new Matrix(L,n,n);
+ }
+
+ /** Solve A*X = B
+ @param B A Matrix with as many rows as A and any number of columns.
+ @return X so that L*L'*X = B
+ @exception IllegalArgumentException Matrix row dimensions must agree.
+ @exception RuntimeException Matrix is not symmetric positive definite.
+ */
+
+ public Matrix solve (Matrix B) {
+ if (B.getRowDimension() != n) {
+ throw new IllegalArgumentException("Matrix row dimensions must agree.");
+ }
+ if (!isspd) {
+ throw new RuntimeException("Matrix is not symmetric positive definite.");
+ }
+
+ // Copy right hand side.
+ double[][] X = B.getArrayCopy();
+ int nx = B.getColumnDimension();
+
+ // Solve L*Y = B;
+ for (int k = 0; k < n; k++) {
+ for (int j = 0; j < nx; j++) {
+ for (int i = 0; i < k ; i++) {
+ X[k][j] -= X[i][j]*L[k][i];
+ }
+ X[k][j] /= L[k][k];
+ }
+ }
+
+ // Solve L'*X = Y;
+ for (int k = n-1; k >= 0; k--) {
+ for (int j = 0; j < nx; j++) {
+ for (int i = k+1; i < n ; i++) {
+ X[k][j] -= X[i][j]*L[i][k];
+ }
+ X[k][j] /= L[k][k];
+ }
+ }
+
+
+ return new Matrix(X,n,nx);
+ }
+ private static final long serialVersionUID = 1;
+
+}
+
diff --git a/Java/src/Jama/EigenvalueDecomposition.java b/Java/src/Jama/EigenvalueDecomposition.java
new file mode 100644
index 00000000..6657c3fa
--- /dev/null
+++ b/Java/src/Jama/EigenvalueDecomposition.java
@@ -0,0 +1,956 @@
+package Jama;
+import Jama.util.*;
+
+/** Eigenvalues and eigenvectors of a real matrix.
+
+ If A is symmetric, then A = V*D*V' where the eigenvalue matrix D is
+ diagonal and the eigenvector matrix V is orthogonal.
+ I.e. A = V.times(D.times(V.transpose())) and
+ V.times(V.transpose()) equals the identity matrix.
+
+ If A is not symmetric, then the eigenvalue matrix D is block diagonal
+ with the real eigenvalues in 1-by-1 blocks and any complex eigenvalues,
+ lambda + i*mu, in 2-by-2 blocks, [lambda, mu; -mu, lambda]. The
+ columns of V represent the eigenvectors in the sense that A*V = V*D,
+ i.e. A.times(V) equals V.times(D). The matrix V may be badly
+ conditioned, or even singular, so the validity of the equation
+ A = V*D*inverse(V) depends upon V.cond().
+**/
+
+public class EigenvalueDecomposition implements java.io.Serializable {
+
+/* ------------------------
+ Class variables
+ * ------------------------ */
+
+ /** Row and column dimension (square matrix).
+ @serial matrix dimension.
+ */
+ private int n;
+
+ /** Symmetry flag.
+ @serial internal symmetry flag.
+ */
+ private boolean issymmetric;
+
+ /** Arrays for internal storage of eigenvalues.
+ @serial internal storage of eigenvalues.
+ */
+ private double[] d, e;
+
+ /** Array for internal storage of eigenvectors.
+ @serial internal storage of eigenvectors.
+ */
+ private double[][] V;
+
+ /** Array for internal storage of nonsymmetric Hessenberg form.
+ @serial internal storage of nonsymmetric Hessenberg form.
+ */
+ private double[][] H;
+
+ /** Working storage for nonsymmetric algorithm.
+ @serial working storage for nonsymmetric algorithm.
+ */
+ private double[] ort;
+
+/* ------------------------
+ Private Methods
+ * ------------------------ */
+
+ // Symmetric Householder reduction to tridiagonal form.
+
+ private void tred2 () {
+
+ // This is derived from the Algol procedures tred2 by
+ // Bowdler, Martin, Reinsch, and Wilkinson, Handbook for
+ // Auto. Comp., Vol.ii-Linear Algebra, and the corresponding
+ // Fortran subroutine in EISPACK.
+
+ for (int j = 0; j < n; j++) {
+ d[j] = V[n-1][j];
+ }
+
+ // Householder reduction to tridiagonal form.
+
+ for (int i = n-1; i > 0; i--) {
+
+ // Scale to avoid under/overflow.
+
+ double scale = 0.0;
+ double h = 0.0;
+ for (int k = 0; k < i; k++) {
+ scale = scale + Math.abs(d[k]);
+ }
+ if (scale == 0.0) {
+ e[i] = d[i-1];
+ for (int j = 0; j < i; j++) {
+ d[j] = V[i-1][j];
+ V[i][j] = 0.0;
+ V[j][i] = 0.0;
+ }
+ } else {
+
+ // Generate Householder vector.
+
+ for (int k = 0; k < i; k++) {
+ d[k] /= scale;
+ h += d[k] * d[k];
+ }
+ double f = d[i-1];
+ double g = Math.sqrt(h);
+ if (f > 0) {
+ g = -g;
+ }
+ e[i] = scale * g;
+ h = h - f * g;
+ d[i-1] = f - g;
+ for (int j = 0; j < i; j++) {
+ e[j] = 0.0;
+ }
+
+ // Apply similarity transformation to remaining columns.
+
+ for (int j = 0; j < i; j++) {
+ f = d[j];
+ V[j][i] = f;
+ g = e[j] + V[j][j] * f;
+ for (int k = j+1; k <= i-1; k++) {
+ g += V[k][j] * d[k];
+ e[k] += V[k][j] * f;
+ }
+ e[j] = g;
+ }
+ f = 0.0;
+ for (int j = 0; j < i; j++) {
+ e[j] /= h;
+ f += e[j] * d[j];
+ }
+ double hh = f / (h + h);
+ for (int j = 0; j < i; j++) {
+ e[j] -= hh * d[j];
+ }
+ for (int j = 0; j < i; j++) {
+ f = d[j];
+ g = e[j];
+ for (int k = j; k <= i-1; k++) {
+ V[k][j] -= (f * e[k] + g * d[k]);
+ }
+ d[j] = V[i-1][j];
+ V[i][j] = 0.0;
+ }
+ }
+ d[i] = h;
+ }
+
+ // Accumulate transformations.
+
+ for (int i = 0; i < n-1; i++) {
+ V[n-1][i] = V[i][i];
+ V[i][i] = 1.0;
+ double h = d[i+1];
+ if (h != 0.0) {
+ for (int k = 0; k <= i; k++) {
+ d[k] = V[k][i+1] / h;
+ }
+ for (int j = 0; j <= i; j++) {
+ double g = 0.0;
+ for (int k = 0; k <= i; k++) {
+ g += V[k][i+1] * V[k][j];
+ }
+ for (int k = 0; k <= i; k++) {
+ V[k][j] -= g * d[k];
+ }
+ }
+ }
+ for (int k = 0; k <= i; k++) {
+ V[k][i+1] = 0.0;
+ }
+ }
+ for (int j = 0; j < n; j++) {
+ d[j] = V[n-1][j];
+ V[n-1][j] = 0.0;
+ }
+ V[n-1][n-1] = 1.0;
+ e[0] = 0.0;
+ }
+
+ // Symmetric tridiagonal QL algorithm.
+
+ private void tql2 () {
+
+ // This is derived from the Algol procedures tql2, by
+ // Bowdler, Martin, Reinsch, and Wilkinson, Handbook for
+ // Auto. Comp., Vol.ii-Linear Algebra, and the corresponding
+ // Fortran subroutine in EISPACK.
+
+ for (int i = 1; i < n; i++) {
+ e[i-1] = e[i];
+ }
+ e[n-1] = 0.0;
+
+ double f = 0.0;
+ double tst1 = 0.0;
+ double eps = Math.pow(2.0,-52.0);
+ for (int l = 0; l < n; l++) {
+
+ // Find small subdiagonal element
+
+ tst1 = Math.max(tst1,Math.abs(d[l]) + Math.abs(e[l]));
+ int m = l;
+ while (m < n) {
+ if (Math.abs(e[m]) <= eps*tst1) {
+ break;
+ }
+ m++;
+ }
+
+ // If m == l, d[l] is an eigenvalue,
+ // otherwise, iterate.
+
+ if (m > l) {
+ int iter = 0;
+ do {
+ iter = iter + 1; // (Could check iteration count here.)
+
+ // Compute implicit shift
+
+ double g = d[l];
+ double p = (d[l+1] - g) / (2.0 * e[l]);
+ double r = Maths.hypot(p,1.0);
+ if (p < 0) {
+ r = -r;
+ }
+ d[l] = e[l] / (p + r);
+ d[l+1] = e[l] * (p + r);
+ double dl1 = d[l+1];
+ double h = g - d[l];
+ for (int i = l+2; i < n; i++) {
+ d[i] -= h;
+ }
+ f = f + h;
+
+ // Implicit QL transformation.
+
+ p = d[m];
+ double c = 1.0;
+ double c2 = c;
+ double c3 = c;
+ double el1 = e[l+1];
+ double s = 0.0;
+ double s2 = 0.0;
+ for (int i = m-1; i >= l; i--) {
+ c3 = c2;
+ c2 = c;
+ s2 = s;
+ g = c * e[i];
+ h = c * p;
+ r = Maths.hypot(p,e[i]);
+ e[i+1] = s * r;
+ s = e[i] / r;
+ c = p / r;
+ p = c * d[i] - s * g;
+ d[i+1] = h + s * (c * g + s * d[i]);
+
+ // Accumulate transformation.
+
+ for (int k = 0; k < n; k++) {
+ h = V[k][i+1];
+ V[k][i+1] = s * V[k][i] + c * h;
+ V[k][i] = c * V[k][i] - s * h;
+ }
+ }
+ p = -s * s2 * c3 * el1 * e[l] / dl1;
+ e[l] = s * p;
+ d[l] = c * p;
+
+ // Check for convergence.
+
+ } while (Math.abs(e[l]) > eps*tst1);
+ }
+ d[l] = d[l] + f;
+ e[l] = 0.0;
+ }
+
+ // Sort eigenvalues and corresponding vectors.
+
+ for (int i = 0; i < n-1; i++) {
+ int k = i;
+ double p = d[i];
+ for (int j = i+1; j < n; j++) {
+ if (d[j] < p) {
+ k = j;
+ p = d[j];
+ }
+ }
+ if (k != i) {
+ d[k] = d[i];
+ d[i] = p;
+ for (int j = 0; j < n; j++) {
+ p = V[j][i];
+ V[j][i] = V[j][k];
+ V[j][k] = p;
+ }
+ }
+ }
+ }
+
+ // Nonsymmetric reduction to Hessenberg form.
+
+ private void orthes () {
+
+ // This is derived from the Algol procedures orthes and ortran,
+ // by Martin and Wilkinson, Handbook for Auto. Comp.,
+ // Vol.ii-Linear Algebra, and the corresponding
+ // Fortran subroutines in EISPACK.
+
+ int low = 0;
+ int high = n-1;
+
+ for (int m = low+1; m <= high-1; m++) {
+
+ // Scale column.
+
+ double scale = 0.0;
+ for (int i = m; i <= high; i++) {
+ scale = scale + Math.abs(H[i][m-1]);
+ }
+ if (scale != 0.0) {
+
+ // Compute Householder transformation.
+
+ double h = 0.0;
+ for (int i = high; i >= m; i--) {
+ ort[i] = H[i][m-1]/scale;
+ h += ort[i] * ort[i];
+ }
+ double g = Math.sqrt(h);
+ if (ort[m] > 0) {
+ g = -g;
+ }
+ h = h - ort[m] * g;
+ ort[m] = ort[m] - g;
+
+ // Apply Householder similarity transformation
+ // H = (I-u*u'/h)*H*(I-u*u')/h)
+
+ for (int j = m; j < n; j++) {
+ double f = 0.0;
+ for (int i = high; i >= m; i--) {
+ f += ort[i]*H[i][j];
+ }
+ f = f/h;
+ for (int i = m; i <= high; i++) {
+ H[i][j] -= f*ort[i];
+ }
+ }
+
+ for (int i = 0; i <= high; i++) {
+ double f = 0.0;
+ for (int j = high; j >= m; j--) {
+ f += ort[j]*H[i][j];
+ }
+ f = f/h;
+ for (int j = m; j <= high; j++) {
+ H[i][j] -= f*ort[j];
+ }
+ }
+ ort[m] = scale*ort[m];
+ H[m][m-1] = scale*g;
+ }
+ }
+
+ // Accumulate transformations (Algol's ortran).
+
+ for (int i = 0; i < n; i++) {
+ for (int j = 0; j < n; j++) {
+ V[i][j] = (i == j ? 1.0 : 0.0);
+ }
+ }
+
+ for (int m = high-1; m >= low+1; m--) {
+ if (H[m][m-1] != 0.0) {
+ for (int i = m+1; i <= high; i++) {
+ ort[i] = H[i][m-1];
+ }
+ for (int j = m; j <= high; j++) {
+ double g = 0.0;
+ for (int i = m; i <= high; i++) {
+ g += ort[i] * V[i][j];
+ }
+ // Double division avoids possible underflow
+ g = (g / ort[m]) / H[m][m-1];
+ for (int i = m; i <= high; i++) {
+ V[i][j] += g * ort[i];
+ }
+ }
+ }
+ }
+ }
+
+
+ // Complex scalar division.
+
+ private transient double cdivr, cdivi;
+ private void cdiv(double xr, double xi, double yr, double yi) {
+ double r,d;
+ if (Math.abs(yr) > Math.abs(yi)) {
+ r = yi/yr;
+ d = yr + r*yi;
+ cdivr = (xr + r*xi)/d;
+ cdivi = (xi - r*xr)/d;
+ } else {
+ r = yr/yi;
+ d = yi + r*yr;
+ cdivr = (r*xr + xi)/d;
+ cdivi = (r*xi - xr)/d;
+ }
+ }
+
+
+ // Nonsymmetric reduction from Hessenberg to real Schur form.
+
+ private void hqr2 () {
+
+ // This is derived from the Algol procedure hqr2,
+ // by Martin and Wilkinson, Handbook for Auto. Comp.,
+ // Vol.ii-Linear Algebra, and the corresponding
+ // Fortran subroutine in EISPACK.
+
+ // Initialize
+
+ int nn = this.n;
+ int n = nn-1;
+ int low = 0;
+ int high = nn-1;
+ double eps = Math.pow(2.0,-52.0);
+ double exshift = 0.0;
+ double p=0,q=0,r=0,s=0,z=0,t,w,x,y;
+
+ // Store roots isolated by balanc and compute matrix norm
+
+ double norm = 0.0;
+ for (int i = 0; i < nn; i++) {
+ if (i < low | i > high) {
+ d[i] = H[i][i];
+ e[i] = 0.0;
+ }
+ for (int j = Math.max(i-1,0); j < nn; j++) {
+ norm = norm + Math.abs(H[i][j]);
+ }
+ }
+
+ // Outer loop over eigenvalue index
+
+ int iter = 0;
+ while (n >= low) {
+
+ // Look for single small sub-diagonal element
+
+ int l = n;
+ while (l > low) {
+ s = Math.abs(H[l-1][l-1]) + Math.abs(H[l][l]);
+ if (s == 0.0) {
+ s = norm;
+ }
+ if (Math.abs(H[l][l-1]) < eps * s) {
+ break;
+ }
+ l--;
+ }
+
+ // Check for convergence
+ // One root found
+
+ if (l == n) {
+ H[n][n] = H[n][n] + exshift;
+ d[n] = H[n][n];
+ e[n] = 0.0;
+ n--;
+ iter = 0;
+
+ // Two roots found
+
+ } else if (l == n-1) {
+ w = H[n][n-1] * H[n-1][n];
+ p = (H[n-1][n-1] - H[n][n]) / 2.0;
+ q = p * p + w;
+ z = Math.sqrt(Math.abs(q));
+ H[n][n] = H[n][n] + exshift;
+ H[n-1][n-1] = H[n-1][n-1] + exshift;
+ x = H[n][n];
+
+ // Real pair
+
+ if (q >= 0) {
+ if (p >= 0) {
+ z = p + z;
+ } else {
+ z = p - z;
+ }
+ d[n-1] = x + z;
+ d[n] = d[n-1];
+ if (z != 0.0) {
+ d[n] = x - w / z;
+ }
+ e[n-1] = 0.0;
+ e[n] = 0.0;
+ x = H[n][n-1];
+ s = Math.abs(x) + Math.abs(z);
+ p = x / s;
+ q = z / s;
+ r = Math.sqrt(p * p+q * q);
+ p = p / r;
+ q = q / r;
+
+ // Row modification
+
+ for (int j = n-1; j < nn; j++) {
+ z = H[n-1][j];
+ H[n-1][j] = q * z + p * H[n][j];
+ H[n][j] = q * H[n][j] - p * z;
+ }
+
+ // Column modification
+
+ for (int i = 0; i <= n; i++) {
+ z = H[i][n-1];
+ H[i][n-1] = q * z + p * H[i][n];
+ H[i][n] = q * H[i][n] - p * z;
+ }
+
+ // Accumulate transformations
+
+ for (int i = low; i <= high; i++) {
+ z = V[i][n-1];
+ V[i][n-1] = q * z + p * V[i][n];
+ V[i][n] = q * V[i][n] - p * z;
+ }
+
+ // Complex pair
+
+ } else {
+ d[n-1] = x + p;
+ d[n] = x + p;
+ e[n-1] = z;
+ e[n] = -z;
+ }
+ n = n - 2;
+ iter = 0;
+
+ // No convergence yet
+
+ } else {
+
+ // Form shift
+
+ x = H[n][n];
+ y = 0.0;
+ w = 0.0;
+ if (l < n) {
+ y = H[n-1][n-1];
+ w = H[n][n-1] * H[n-1][n];
+ }
+
+ // Wilkinson's original ad hoc shift
+
+ if (iter == 10) {
+ exshift += x;
+ for (int i = low; i <= n; i++) {
+ H[i][i] -= x;
+ }
+ s = Math.abs(H[n][n-1]) + Math.abs(H[n-1][n-2]);
+ x = y = 0.75 * s;
+ w = -0.4375 * s * s;
+ }
+
+ // MATLAB's new ad hoc shift
+
+ if (iter == 30) {
+ s = (y - x) / 2.0;
+ s = s * s + w;
+ if (s > 0) {
+ s = Math.sqrt(s);
+ if (y < x) {
+ s = -s;
+ }
+ s = x - w / ((y - x) / 2.0 + s);
+ for (int i = low; i <= n; i++) {
+ H[i][i] -= s;
+ }
+ exshift += s;
+ x = y = w = 0.964;
+ }
+ }
+
+ iter = iter + 1; // (Could check iteration count here.)
+
+ // Look for two consecutive small sub-diagonal elements
+
+ int m = n-2;
+ while (m >= l) {
+ z = H[m][m];
+ r = x - z;
+ s = y - z;
+ p = (r * s - w) / H[m+1][m] + H[m][m+1];
+ q = H[m+1][m+1] - z - r - s;
+ r = H[m+2][m+1];
+ s = Math.abs(p) + Math.abs(q) + Math.abs(r);
+ p = p / s;
+ q = q / s;
+ r = r / s;
+ if (m == l) {
+ break;
+ }
+ if (Math.abs(H[m][m-1]) * (Math.abs(q) + Math.abs(r)) <
+ eps * (Math.abs(p) * (Math.abs(H[m-1][m-1]) + Math.abs(z) +
+ Math.abs(H[m+1][m+1])))) {
+ break;
+ }
+ m--;
+ }
+
+ for (int i = m+2; i <= n; i++) {
+ H[i][i-2] = 0.0;
+ if (i > m+2) {
+ H[i][i-3] = 0.0;
+ }
+ }
+
+ // Double QR step involving rows l:n and columns m:n
+
+
+ for (int k = m; k <= n-1; k++) {
+ boolean notlast = (k != n-1);
+ if (k != m) {
+ p = H[k][k-1];
+ q = H[k+1][k-1];
+ r = (notlast ? H[k+2][k-1] : 0.0);
+ x = Math.abs(p) + Math.abs(q) + Math.abs(r);
+ if (x == 0.0) {
+ continue;
+ }
+ p = p / x;
+ q = q / x;
+ r = r / x;
+ }
+
+ s = Math.sqrt(p * p + q * q + r * r);
+ if (p < 0) {
+ s = -s;
+ }
+ if (s != 0) {
+ if (k != m) {
+ H[k][k-1] = -s * x;
+ } else if (l != m) {
+ H[k][k-1] = -H[k][k-1];
+ }
+ p = p + s;
+ x = p / s;
+ y = q / s;
+ z = r / s;
+ q = q / p;
+ r = r / p;
+
+ // Row modification
+
+ for (int j = k; j < nn; j++) {
+ p = H[k][j] + q * H[k+1][j];
+ if (notlast) {
+ p = p + r * H[k+2][j];
+ H[k+2][j] = H[k+2][j] - p * z;
+ }
+ H[k][j] = H[k][j] - p * x;
+ H[k+1][j] = H[k+1][j] - p * y;
+ }
+
+ // Column modification
+
+ for (int i = 0; i <= Math.min(n,k+3); i++) {
+ p = x * H[i][k] + y * H[i][k+1];
+ if (notlast) {
+ p = p + z * H[i][k+2];
+ H[i][k+2] = H[i][k+2] - p * r;
+ }
+ H[i][k] = H[i][k] - p;
+ H[i][k+1] = H[i][k+1] - p * q;
+ }
+
+ // Accumulate transformations
+
+ for (int i = low; i <= high; i++) {
+ p = x * V[i][k] + y * V[i][k+1];
+ if (notlast) {
+ p = p + z * V[i][k+2];
+ V[i][k+2] = V[i][k+2] - p * r;
+ }
+ V[i][k] = V[i][k] - p;
+ V[i][k+1] = V[i][k+1] - p * q;
+ }
+ } // (s != 0)
+ } // k loop
+ } // check convergence
+ } // while (n >= low)
+
+ // Backsubstitute to find vectors of upper triangular form
+
+ if (norm == 0.0) {
+ return;
+ }
+
+ for (n = nn-1; n >= 0; n--) {
+ p = d[n];
+ q = e[n];
+
+ // Real vector
+
+ if (q == 0) {
+ int l = n;
+ H[n][n] = 1.0;
+ for (int i = n-1; i >= 0; i--) {
+ w = H[i][i] - p;
+ r = 0.0;
+ for (int j = l; j <= n; j++) {
+ r = r + H[i][j] * H[j][n];
+ }
+ if (e[i] < 0.0) {
+ z = w;
+ s = r;
+ } else {
+ l = i;
+ if (e[i] == 0.0) {
+ if (w != 0.0) {
+ H[i][n] = -r / w;
+ } else {
+ H[i][n] = -r / (eps * norm);
+ }
+
+ // Solve real equations
+
+ } else {
+ x = H[i][i+1];
+ y = H[i+1][i];
+ q = (d[i] - p) * (d[i] - p) + e[i] * e[i];
+ t = (x * s - z * r) / q;
+ H[i][n] = t;
+ if (Math.abs(x) > Math.abs(z)) {
+ H[i+1][n] = (-r - w * t) / x;
+ } else {
+ H[i+1][n] = (-s - y * t) / z;
+ }
+ }
+
+ // Overflow control
+
+ t = Math.abs(H[i][n]);
+ if ((eps * t) * t > 1) {
+ for (int j = i; j <= n; j++) {
+ H[j][n] = H[j][n] / t;
+ }
+ }
+ }
+ }
+
+ // Complex vector
+
+ } else if (q < 0) {
+ int l = n-1;
+
+ // Last vector component imaginary so matrix is triangular
+
+ if (Math.abs(H[n][n-1]) > Math.abs(H[n-1][n])) {
+ H[n-1][n-1] = q / H[n][n-1];
+ H[n-1][n] = -(H[n][n] - p) / H[n][n-1];
+ } else {
+ cdiv(0.0,-H[n-1][n],H[n-1][n-1]-p,q);
+ H[n-1][n-1] = cdivr;
+ H[n-1][n] = cdivi;
+ }
+ H[n][n-1] = 0.0;
+ H[n][n] = 1.0;
+ for (int i = n-2; i >= 0; i--) {
+ double ra,sa,vr,vi;
+ ra = 0.0;
+ sa = 0.0;
+ for (int j = l; j <= n; j++) {
+ ra = ra + H[i][j] * H[j][n-1];
+ sa = sa + H[i][j] * H[j][n];
+ }
+ w = H[i][i] - p;
+
+ if (e[i] < 0.0) {
+ z = w;
+ r = ra;
+ s = sa;
+ } else {
+ l = i;
+ if (e[i] == 0) {
+ cdiv(-ra,-sa,w,q);
+ H[i][n-1] = cdivr;
+ H[i][n] = cdivi;
+ } else {
+
+ // Solve complex equations
+
+ x = H[i][i+1];
+ y = H[i+1][i];
+ vr = (d[i] - p) * (d[i] - p) + e[i] * e[i] - q * q;
+ vi = (d[i] - p) * 2.0 * q;
+ if (vr == 0.0 & vi == 0.0) {
+ vr = eps * norm * (Math.abs(w) + Math.abs(q) +
+ Math.abs(x) + Math.abs(y) + Math.abs(z));
+ }
+ cdiv(x*r-z*ra+q*sa,x*s-z*sa-q*ra,vr,vi);
+ H[i][n-1] = cdivr;
+ H[i][n] = cdivi;
+ if (Math.abs(x) > (Math.abs(z) + Math.abs(q))) {
+ H[i+1][n-1] = (-ra - w * H[i][n-1] + q * H[i][n]) / x;
+ H[i+1][n] = (-sa - w * H[i][n] - q * H[i][n-1]) / x;
+ } else {
+ cdiv(-r-y*H[i][n-1],-s-y*H[i][n],z,q);
+ H[i+1][n-1] = cdivr;
+ H[i+1][n] = cdivi;
+ }
+ }
+
+ // Overflow control
+
+ t = Math.max(Math.abs(H[i][n-1]),Math.abs(H[i][n]));
+ if ((eps * t) * t > 1) {
+ for (int j = i; j <= n; j++) {
+ H[j][n-1] = H[j][n-1] / t;
+ H[j][n] = H[j][n] / t;
+ }
+ }
+ }
+ }
+ }
+ }
+
+ // Vectors of isolated roots
+
+ for (int i = 0; i < nn; i++) {
+ if (i < low | i > high) {
+ for (int j = i; j < nn; j++) {
+ V[i][j] = H[i][j];
+ }
+ }
+ }
+
+ // Back transformation to get eigenvectors of original matrix
+
+ for (int j = nn-1; j >= low; j--) {
+ for (int i = low; i <= high; i++) {
+ z = 0.0;
+ for (int k = low; k <= Math.min(j,high); k++) {
+ z = z + V[i][k] * H[k][j];
+ }
+ V[i][j] = z;
+ }
+ }
+ }
+
+
+/* ------------------------
+ Constructor
+ * ------------------------ */
+
+ /** Check for symmetry, then construct the eigenvalue decomposition
+ Structure to access D and V.
+ @param Arg Square matrix
+ */
+
+ public EigenvalueDecomposition (Matrix Arg) {
+ double[][] A = Arg.getArray();
+ n = Arg.getColumnDimension();
+ V = new double[n][n];
+ d = new double[n];
+ e = new double[n];
+
+ issymmetric = true;
+ for (int j = 0; (j < n) & issymmetric; j++) {
+ for (int i = 0; (i < n) & issymmetric; i++) {
+ issymmetric = (A[i][j] == A[j][i]);
+ }
+ }
+
+ if (issymmetric) {
+ for (int i = 0; i < n; i++) {
+ for (int j = 0; j < n; j++) {
+ V[i][j] = A[i][j];
+ }
+ }
+
+ // Tridiagonalize.
+ tred2();
+
+ // Diagonalize.
+ tql2();
+
+ } else {
+ H = new double[n][n];
+ ort = new double[n];
+
+ for (int j = 0; j < n; j++) {
+ for (int i = 0; i < n; i++) {
+ H[i][j] = A[i][j];
+ }
+ }
+
+ // Reduce to Hessenberg form.
+ orthes();
+
+ // Reduce Hessenberg to real Schur form.
+ hqr2();
+ }
+ }
+
+/* ------------------------
+ Public Methods
+ * ------------------------ */
+
+ /** Return the eigenvector matrix
+ @return V
+ */
+
+ public Matrix getV () {
+ return new Matrix(V,n,n);
+ }
+
+ /** Return the real parts of the eigenvalues
+ @return real(diag(D))
+ */
+
+ public double[] getRealEigenvalues () {
+ return d;
+ }
+
+ /** Return the imaginary parts of the eigenvalues
+ @return imag(diag(D))
+ */
+
+ public double[] getImagEigenvalues () {
+ return e;
+ }
+
+ /** Return the block diagonal eigenvalue matrix
+ @return D
+ */
+
+ public Matrix getD () {
+ Matrix X = new Matrix(n,n);
+ double[][] D = X.getArray();
+ for (int i = 0; i < n; i++) {
+ for (int j = 0; j < n; j++) {
+ D[i][j] = 0.0;
+ }
+ D[i][i] = d[i];
+ if (e[i] > 0) {
+ D[i][i+1] = e[i];
+ } else if (e[i] < 0) {
+ D[i][i-1] = e[i];
+ }
+ }
+ return X;
+ }
+ private static final long serialVersionUID = 1;
+}
diff --git a/Java/src/Jama/LUDecomposition.java b/Java/src/Jama/LUDecomposition.java
new file mode 100644
index 00000000..dc43dc78
--- /dev/null
+++ b/Java/src/Jama/LUDecomposition.java
@@ -0,0 +1,312 @@
+package Jama;
+
+ /** LU Decomposition.
+
+ For an m-by-n matrix A with m >= n, the LU decomposition is an m-by-n
+ unit lower triangular matrix L, an n-by-n upper triangular matrix U,
+ and a permutation vector piv of length m so that A(piv,:) = L*U.
+ If m < n, then L is m-by-m and U is m-by-n.
+
+ The LU decompostion with pivoting always exists, even if the matrix is
+ singular, so the constructor will never fail. The primary use of the
+ LU decomposition is in the solution of square systems of simultaneous
+ linear equations. This will fail if isNonsingular() returns false.
+ */
+
+public class LUDecomposition implements java.io.Serializable {
+
+/* ------------------------
+ Class variables
+ * ------------------------ */
+
+ /** Array for internal storage of decomposition.
+ @serial internal array storage.
+ */
+ private double[][] LU;
+
+ /** Row and column dimensions, and pivot sign.
+ @serial column dimension.
+ @serial row dimension.
+ @serial pivot sign.
+ */
+ private int m, n, pivsign;
+
+ /** Internal storage of pivot vector.
+ @serial pivot vector.
+ */
+ private int[] piv;
+
+/* ------------------------
+ Constructor
+ * ------------------------ */
+
+ /** LU Decomposition
+ Structure to access L, U and piv.
+ @param A Rectangular matrix
+ */
+
+ public LUDecomposition (Matrix A) {
+
+ // Use a "left-looking", dot-product, Crout/Doolittle algorithm.
+
+ LU = A.getArrayCopy();
+ m = A.getRowDimension();
+ n = A.getColumnDimension();
+ piv = new int[m];
+ for (int i = 0; i < m; i++) {
+ piv[i] = i;
+ }
+ pivsign = 1;
+ double[] LUrowi;
+ double[] LUcolj = new double[m];
+
+ // Outer loop.
+
+ for (int j = 0; j < n; j++) {
+
+ // Make a copy of the j-th column to localize references.
+
+ for (int i = 0; i < m; i++) {
+ LUcolj[i] = LU[i][j];
+ }
+
+ // Apply previous transformations.
+
+ for (int i = 0; i < m; i++) {
+ LUrowi = LU[i];
+
+ // Most of the time is spent in the following dot product.
+
+ int kmax = Math.min(i,j);
+ double s = 0.0;
+ for (int k = 0; k < kmax; k++) {
+ s += LUrowi[k]*LUcolj[k];
+ }
+
+ LUrowi[j] = LUcolj[i] -= s;
+ }
+
+ // Find pivot and exchange if necessary.
+
+ int p = j;
+ for (int i = j+1; i < m; i++) {
+ if (Math.abs(LUcolj[i]) > Math.abs(LUcolj[p])) {
+ p = i;
+ }
+ }
+ if (p != j) {
+ for (int k = 0; k < n; k++) {
+ double t = LU[p][k]; LU[p][k] = LU[j][k]; LU[j][k] = t;
+ }
+ int k = piv[p]; piv[p] = piv[j]; piv[j] = k;
+ pivsign = -pivsign;
+ }
+
+ // Compute multipliers.
+
+ if (j < m & LU[j][j] != 0.0) {
+ for (int i = j+1; i < m; i++) {
+ LU[i][j] /= LU[j][j];
+ }
+ }
+ }
+ }
+
+/* ------------------------
+ Temporary, experimental code.
+ ------------------------ *\
+
+ \** LU Decomposition, computed by Gaussian elimination.
+
+ This constructor computes L and U with the "daxpy"-based elimination
+ algorithm used in LINPACK and MATLAB. In Java, we suspect the dot-product,
+ Crout algorithm will be faster. We have temporarily included this
+ constructor until timing experiments confirm this suspicion.
+
+ @param A Rectangular matrix
+ @param linpackflag Use Gaussian elimination. Actual value ignored.
+ @return Structure to access L, U and piv.
+ *\
+
+ public LUDecomposition (Matrix A, int linpackflag) {
+ // Initialize.
+ LU = A.getArrayCopy();
+ m = A.getRowDimension();
+ n = A.getColumnDimension();
+ piv = new int[m];
+ for (int i = 0; i < m; i++) {
+ piv[i] = i;
+ }
+ pivsign = 1;
+ // Main loop.
+ for (int k = 0; k < n; k++) {
+ // Find pivot.
+ int p = k;
+ for (int i = k+1; i < m; i++) {
+ if (Math.abs(LU[i][k]) > Math.abs(LU[p][k])) {
+ p = i;
+ }
+ }
+ // Exchange if necessary.
+ if (p != k) {
+ for (int j = 0; j < n; j++) {
+ double t = LU[p][j]; LU[p][j] = LU[k][j]; LU[k][j] = t;
+ }
+ int t = piv[p]; piv[p] = piv[k]; piv[k] = t;
+ pivsign = -pivsign;
+ }
+ // Compute multipliers and eliminate k-th column.
+ if (LU[k][k] != 0.0) {
+ for (int i = k+1; i < m; i++) {
+ LU[i][k] /= LU[k][k];
+ for (int j = k+1; j < n; j++) {
+ LU[i][j] -= LU[i][k]*LU[k][j];
+ }
+ }
+ }
+ }
+ }
+
+\* ------------------------
+ End of temporary code.
+ * ------------------------ */
+
+/* ------------------------
+ Public Methods
+ * ------------------------ */
+
+ /** Is the matrix nonsingular?
+ @return true if U, and hence A, is nonsingular.
+ */
+
+ public boolean isNonsingular () {
+ for (int j = 0; j < n; j++) {
+ if (LU[j][j] == 0)
+ return false;
+ }
+ return true;
+ }
+
+ /** Return lower triangular factor
+ @return L
+ */
+
+ public Matrix getL () {
+ Matrix X = new Matrix(m,n);
+ double[][] L = X.getArray();
+ for (int i = 0; i < m; i++) {
+ for (int j = 0; j < n; j++) {
+ if (i > j) {
+ L[i][j] = LU[i][j];
+ } else if (i == j) {
+ L[i][j] = 1.0;
+ } else {
+ L[i][j] = 0.0;
+ }
+ }
+ }
+ return X;
+ }
+
+ /** Return upper triangular factor
+ @return U
+ */
+
+ public Matrix getU () {
+ Matrix X = new Matrix(n,n);
+ double[][] U = X.getArray();
+ for (int i = 0; i < n; i++) {
+ for (int j = 0; j < n; j++) {
+ if (i <= j) {
+ U[i][j] = LU[i][j];
+ } else {
+ U[i][j] = 0.0;
+ }
+ }
+ }
+ return X;
+ }
+
+ /** Return pivot permutation vector
+ @return piv
+ */
+
+ public int[] getPivot () {
+ int[] p = new int[m];
+ for (int i = 0; i < m; i++) {
+ p[i] = piv[i];
+ }
+ return p;
+ }
+
+ /** Return pivot permutation vector as a one-dimensional double array
+ @return (double) piv
+ */
+
+ public double[] getDoublePivot () {
+ double[] vals = new double[m];
+ for (int i = 0; i < m; i++) {
+ vals[i] = (double) piv[i];
+ }
+ return vals;
+ }
+
+ /** Determinant
+ @return det(A)
+ @exception IllegalArgumentException Matrix must be square
+ */
+
+ public double det () {
+ if (m != n) {
+ throw new IllegalArgumentException("Matrix must be square.");
+ }
+ double d = (double) pivsign;
+ for (int j = 0; j < n; j++) {
+ d *= LU[j][j];
+ }
+ return d;
+ }
+
+ /** Solve A*X = B
+ @param B A Matrix with as many rows as A and any number of columns.
+ @return X so that L*U*X = B(piv,:)
+ @exception IllegalArgumentException Matrix row dimensions must agree.
+ @exception RuntimeException Matrix is singular.
+ */
+
+ public Matrix solve (Matrix B) {
+ if (B.getRowDimension() != m) {
+ throw new IllegalArgumentException("Matrix row dimensions must agree.");
+ }
+ if (!this.isNonsingular()) {
+ throw new RuntimeException("Matrix is singular.");
+ }
+
+ // Copy right hand side with pivoting
+ int nx = B.getColumnDimension();
+ Matrix Xmat = B.getMatrix(piv,0,nx-1);
+ double[][] X = Xmat.getArray();
+
+ // Solve L*Y = B(piv,:)
+ for (int k = 0; k < n; k++) {
+ for (int i = k+1; i < n; i++) {
+ for (int j = 0; j < nx; j++) {
+ X[i][j] -= X[k][j]*LU[i][k];
+ }
+ }
+ }
+ // Solve U*X = Y;
+ for (int k = n-1; k >= 0; k--) {
+ for (int j = 0; j < nx; j++) {
+ X[k][j] /= LU[k][k];
+ }
+ for (int i = 0; i < k; i++) {
+ for (int j = 0; j < nx; j++) {
+ X[i][j] -= X[k][j]*LU[i][k];
+ }
+ }
+ }
+ return Xmat;
+ }
+ private static final long serialVersionUID = 1;
+}
diff --git a/Java/src/Jama/Matrix.java b/Java/src/Jama/Matrix.java
new file mode 100644
index 00000000..8b113877
--- /dev/null
+++ b/Java/src/Jama/Matrix.java
@@ -0,0 +1,1046 @@
+package Jama;
+
+import java.text.NumberFormat;
+import java.text.DecimalFormat;
+import java.text.DecimalFormatSymbols;
+import java.util.Locale;
+import java.text.FieldPosition;
+import java.io.PrintWriter;
+import java.io.BufferedReader;
+import java.io.StreamTokenizer;
+import Jama.util.*;
+
+/**
+ Jama = Java Matrix class.
+
+ The Java Matrix Class provides the fundamental operations of numerical
+ linear algebra. Various constructors create Matrices from two dimensional
+ arrays of double precision floating point numbers. Various "gets" and
+ "sets" provide access to submatrices and matrix elements. Several methods
+ implement basic matrix arithmetic, including matrix addition and
+ multiplication, matrix norms, and element-by-element array operations.
+ Methods for reading and printing matrices are also included. All the
+ operations in this version of the Matrix Class involve real matrices.
+ Complex matrices may be handled in a future version.
+
+ Five fundamental matrix decompositions, which consist of pairs or triples
+ of matrices, permutation vectors, and the like, produce results in five
+ decomposition classes. These decompositions are accessed by the Matrix
+ class to compute solutions of simultaneous linear equations, determinants,
+ inverses and other matrix functions. The five decompositions are:
+
+
Cholesky Decomposition of symmetric, positive definite matrices.
+
LU Decomposition of rectangular matrices.
+
QR Decomposition of rectangular matrices.
+
Singular Value Decomposition of rectangular matrices.
+
Eigenvalue Decomposition of both symmetric and nonsymmetric square matrices.
+
+
+
Example of use:
+
+
Solve a linear system A x = b and compute the residual norm, ||b - A x||.
+
+ double[][] vals = {{1.,2.,3},{4.,5.,6.},{7.,8.,10.}};
+ Matrix A = new Matrix(vals);
+ Matrix b = Matrix.random(3,1);
+ Matrix x = A.solve(b);
+ Matrix r = A.times(x).minus(b);
+ double rnorm = r.normInf();
+
+
+
+@author The MathWorks, Inc. and the National Institute of Standards and Technology.
+@version 5 August 1998
+*/
+
+public class Matrix implements Cloneable, java.io.Serializable {
+
+/* ------------------------
+ Class variables
+ * ------------------------ */
+
+ /** Array for internal storage of elements.
+ @serial internal array storage.
+ */
+ private double[][] A;
+
+ /** Row and column dimensions.
+ @serial row dimension.
+ @serial column dimension.
+ */
+ private int m, n;
+
+/* ------------------------
+ Constructors
+ * ------------------------ */
+
+ /** Construct an m-by-n matrix of zeros.
+ @param m Number of rows.
+ @param n Number of colums.
+ */
+
+ public Matrix (int m, int n) {
+ this.m = m;
+ this.n = n;
+ A = new double[m][n];
+ }
+
+ /** Construct an m-by-n constant matrix.
+ @param m Number of rows.
+ @param n Number of colums.
+ @param s Fill the matrix with this scalar value.
+ */
+
+ public Matrix (int m, int n, double s) {
+ this.m = m;
+ this.n = n;
+ A = new double[m][n];
+ for (int i = 0; i < m; i++) {
+ for (int j = 0; j < n; j++) {
+ A[i][j] = s;
+ }
+ }
+ }
+
+ /** Construct a matrix from a 2-D array.
+ @param A Two-dimensional array of doubles.
+ @exception IllegalArgumentException All rows must have the same length
+ @see #constructWithCopy
+ */
+
+ public Matrix (double[][] A) {
+ m = A.length;
+ n = A[0].length;
+ for (int i = 0; i < m; i++) {
+ if (A[i].length != n) {
+ throw new IllegalArgumentException("All rows must have the same length.");
+ }
+ }
+ this.A = A;
+ }
+
+ /** Construct a matrix quickly without checking arguments.
+ @param A Two-dimensional array of doubles.
+ @param m Number of rows.
+ @param n Number of colums.
+ */
+
+ public Matrix (double[][] A, int m, int n) {
+ this.A = A;
+ this.m = m;
+ this.n = n;
+ }
+
+ /** Construct a matrix from a one-dimensional packed array
+ @param vals One-dimensional array of doubles, packed by columns (ala Fortran).
+ @param m Number of rows.
+ @exception IllegalArgumentException Array length must be a multiple of m.
+ */
+
+ public Matrix (double vals[], int m) {
+ this.m = m;
+ n = (m != 0 ? vals.length/m : 0);
+ if (m*n != vals.length) {
+ throw new IllegalArgumentException("Array length must be a multiple of m.");
+ }
+ A = new double[m][n];
+ for (int i = 0; i < m; i++) {
+ for (int j = 0; j < n; j++) {
+ A[i][j] = vals[i+j*m];
+ }
+ }
+ }
+
+/* ------------------------
+ Public Methods
+ * ------------------------ */
+
+ /** Construct a matrix from a copy of a 2-D array.
+ @param A Two-dimensional array of doubles.
+ @exception IllegalArgumentException All rows must have the same length
+ */
+
+ public static Matrix constructWithCopy(double[][] A) {
+ int m = A.length;
+ int n = A[0].length;
+ Matrix X = new Matrix(m,n);
+ double[][] C = X.getArray();
+ for (int i = 0; i < m; i++) {
+ if (A[i].length != n) {
+ throw new IllegalArgumentException
+ ("All rows must have the same length.");
+ }
+ for (int j = 0; j < n; j++) {
+ C[i][j] = A[i][j];
+ }
+ }
+ return X;
+ }
+
+ /** Make a deep copy of a matrix
+ */
+
+ public Matrix copy () {
+ Matrix X = new Matrix(m,n);
+ double[][] C = X.getArray();
+ for (int i = 0; i < m; i++) {
+ for (int j = 0; j < n; j++) {
+ C[i][j] = A[i][j];
+ }
+ }
+ return X;
+ }
+
+ /** Clone the Matrix object.
+ */
+
+ public Object clone () {
+ return this.copy();
+ }
+
+ /** Access the internal two-dimensional array.
+ @return Pointer to the two-dimensional array of matrix elements.
+ */
+
+ public double[][] getArray () {
+ return A;
+ }
+
+ /** Copy the internal two-dimensional array.
+ @return Two-dimensional array copy of matrix elements.
+ */
+
+ public double[][] getArrayCopy () {
+ double[][] C = new double[m][n];
+ for (int i = 0; i < m; i++) {
+ for (int j = 0; j < n; j++) {
+ C[i][j] = A[i][j];
+ }
+ }
+ return C;
+ }
+
+ /** Make a one-dimensional column packed copy of the internal array.
+ @return Matrix elements packed in a one-dimensional array by columns.
+ */
+
+ public double[] getColumnPackedCopy () {
+ double[] vals = new double[m*n];
+ for (int i = 0; i < m; i++) {
+ for (int j = 0; j < n; j++) {
+ vals[i+j*m] = A[i][j];
+ }
+ }
+ return vals;
+ }
+
+ /** Make a one-dimensional row packed copy of the internal array.
+ @return Matrix elements packed in a one-dimensional array by rows.
+ */
+
+ public double[] getRowPackedCopy () {
+ double[] vals = new double[m*n];
+ for (int i = 0; i < m; i++) {
+ for (int j = 0; j < n; j++) {
+ vals[i*n+j] = A[i][j];
+ }
+ }
+ return vals;
+ }
+
+ /** Get row dimension.
+ @return m, the number of rows.
+ */
+
+ public int getRowDimension () {
+ return m;
+ }
+
+ /** Get column dimension.
+ @return n, the number of columns.
+ */
+
+ public int getColumnDimension () {
+ return n;
+ }
+
+ /** Get a single element.
+ @param i Row index.
+ @param j Column index.
+ @return A(i,j)
+ @exception ArrayIndexOutOfBoundsException
+ */
+
+ public double get (int i, int j) {
+ return A[i][j];
+ }
+
+ /** Get a submatrix.
+ @param i0 Initial row index
+ @param i1 Final row index
+ @param j0 Initial column index
+ @param j1 Final column index
+ @return A(i0:i1,j0:j1)
+ @exception ArrayIndexOutOfBoundsException Submatrix indices
+ */
+
+ public Matrix getMatrix (int i0, int i1, int j0, int j1) {
+ Matrix X = new Matrix(i1-i0+1,j1-j0+1);
+ double[][] B = X.getArray();
+ try {
+ for (int i = i0; i <= i1; i++) {
+ for (int j = j0; j <= j1; j++) {
+ B[i-i0][j-j0] = A[i][j];
+ }
+ }
+ } catch(ArrayIndexOutOfBoundsException e) {
+ throw new ArrayIndexOutOfBoundsException("Submatrix indices");
+ }
+ return X;
+ }
+
+ /** Get a submatrix.
+ @param r Array of row indices.
+ @param c Array of column indices.
+ @return A(r(:),c(:))
+ @exception ArrayIndexOutOfBoundsException Submatrix indices
+ */
+
+ public Matrix getMatrix (int[] r, int[] c) {
+ Matrix X = new Matrix(r.length,c.length);
+ double[][] B = X.getArray();
+ try {
+ for (int i = 0; i < r.length; i++) {
+ for (int j = 0; j < c.length; j++) {
+ B[i][j] = A[r[i]][c[j]];
+ }
+ }
+ } catch(ArrayIndexOutOfBoundsException e) {
+ throw new ArrayIndexOutOfBoundsException("Submatrix indices");
+ }
+ return X;
+ }
+
+ /** Get a submatrix.
+ @param i0 Initial row index
+ @param i1 Final row index
+ @param c Array of column indices.
+ @return A(i0:i1,c(:))
+ @exception ArrayIndexOutOfBoundsException Submatrix indices
+ */
+
+ public Matrix getMatrix (int i0, int i1, int[] c) {
+ Matrix X = new Matrix(i1-i0+1,c.length);
+ double[][] B = X.getArray();
+ try {
+ for (int i = i0; i <= i1; i++) {
+ for (int j = 0; j < c.length; j++) {
+ B[i-i0][j] = A[i][c[j]];
+ }
+ }
+ } catch(ArrayIndexOutOfBoundsException e) {
+ throw new ArrayIndexOutOfBoundsException("Submatrix indices");
+ }
+ return X;
+ }
+
+ /** Get a submatrix.
+ @param r Array of row indices.
+ @param j0 Initial column index
+ @param j1 Final column index
+ @return A(r(:),j0:j1)
+ @exception ArrayIndexOutOfBoundsException Submatrix indices
+ */
+
+ public Matrix getMatrix (int[] r, int j0, int j1) {
+ Matrix X = new Matrix(r.length,j1-j0+1);
+ double[][] B = X.getArray();
+ try {
+ for (int i = 0; i < r.length; i++) {
+ for (int j = j0; j <= j1; j++) {
+ B[i][j-j0] = A[r[i]][j];
+ }
+ }
+ } catch(ArrayIndexOutOfBoundsException e) {
+ throw new ArrayIndexOutOfBoundsException("Submatrix indices");
+ }
+ return X;
+ }
+
+ /** Set a single element.
+ @param i Row index.
+ @param j Column index.
+ @param s A(i,j).
+ @exception ArrayIndexOutOfBoundsException
+ */
+
+ public void set (int i, int j, double s) {
+ A[i][j] = s;
+ }
+
+ /** Set a submatrix.
+ @param i0 Initial row index
+ @param i1 Final row index
+ @param j0 Initial column index
+ @param j1 Final column index
+ @param X A(i0:i1,j0:j1)
+ @exception ArrayIndexOutOfBoundsException Submatrix indices
+ */
+
+ public void setMatrix (int i0, int i1, int j0, int j1, Matrix X) {
+ try {
+ for (int i = i0; i <= i1; i++) {
+ for (int j = j0; j <= j1; j++) {
+ A[i][j] = X.get(i-i0,j-j0);
+ }
+ }
+ } catch(ArrayIndexOutOfBoundsException e) {
+ throw new ArrayIndexOutOfBoundsException("Submatrix indices");
+ }
+ }
+
+ /** Set a submatrix.
+ @param r Array of row indices.
+ @param c Array of column indices.
+ @param X A(r(:),c(:))
+ @exception ArrayIndexOutOfBoundsException Submatrix indices
+ */
+
+ public void setMatrix (int[] r, int[] c, Matrix X) {
+ try {
+ for (int i = 0; i < r.length; i++) {
+ for (int j = 0; j < c.length; j++) {
+ A[r[i]][c[j]] = X.get(i,j);
+ }
+ }
+ } catch(ArrayIndexOutOfBoundsException e) {
+ throw new ArrayIndexOutOfBoundsException("Submatrix indices");
+ }
+ }
+
+ /** Set a submatrix.
+ @param r Array of row indices.
+ @param j0 Initial column index
+ @param j1 Final column index
+ @param X A(r(:),j0:j1)
+ @exception ArrayIndexOutOfBoundsException Submatrix indices
+ */
+
+ public void setMatrix (int[] r, int j0, int j1, Matrix X) {
+ try {
+ for (int i = 0; i < r.length; i++) {
+ for (int j = j0; j <= j1; j++) {
+ A[r[i]][j] = X.get(i,j-j0);
+ }
+ }
+ } catch(ArrayIndexOutOfBoundsException e) {
+ throw new ArrayIndexOutOfBoundsException("Submatrix indices");
+ }
+ }
+
+ /** Set a submatrix.
+ @param i0 Initial row index
+ @param i1 Final row index
+ @param c Array of column indices.
+ @param X A(i0:i1,c(:))
+ @exception ArrayIndexOutOfBoundsException Submatrix indices
+ */
+
+ public void setMatrix (int i0, int i1, int[] c, Matrix X) {
+ try {
+ for (int i = i0; i <= i1; i++) {
+ for (int j = 0; j < c.length; j++) {
+ A[i][c[j]] = X.get(i-i0,j);
+ }
+ }
+ } catch(ArrayIndexOutOfBoundsException e) {
+ throw new ArrayIndexOutOfBoundsException("Submatrix indices");
+ }
+ }
+
+ /** Matrix transpose.
+ @return A'
+ */
+
+ public Matrix transpose () {
+ Matrix X = new Matrix(n,m);
+ double[][] C = X.getArray();
+ for (int i = 0; i < m; i++) {
+ for (int j = 0; j < n; j++) {
+ C[j][i] = A[i][j];
+ }
+ }
+ return X;
+ }
+
+ /** One norm
+ @return maximum column sum.
+ */
+
+ public double norm1 () {
+ double f = 0;
+ for (int j = 0; j < n; j++) {
+ double s = 0;
+ for (int i = 0; i < m; i++) {
+ s += Math.abs(A[i][j]);
+ }
+ f = Math.max(f,s);
+ }
+ return f;
+ }
+
+ /** Two norm
+ @return maximum singular value.
+ */
+
+ public double norm2 () {
+ return (new SingularValueDecomposition(this).norm2());
+ }
+
+ /** Infinity norm
+ @return maximum row sum.
+ */
+
+ public double normInf () {
+ double f = 0;
+ for (int i = 0; i < m; i++) {
+ double s = 0;
+ for (int j = 0; j < n; j++) {
+ s += Math.abs(A[i][j]);
+ }
+ f = Math.max(f,s);
+ }
+ return f;
+ }
+
+ /** Frobenius norm
+ @return sqrt of sum of squares of all elements.
+ */
+
+ public double normF () {
+ double f = 0;
+ for (int i = 0; i < m; i++) {
+ for (int j = 0; j < n; j++) {
+ f = Maths.hypot(f,A[i][j]);
+ }
+ }
+ return f;
+ }
+
+ /** Unary minus
+ @return -A
+ */
+
+ public Matrix uminus () {
+ Matrix X = new Matrix(m,n);
+ double[][] C = X.getArray();
+ for (int i = 0; i < m; i++) {
+ for (int j = 0; j < n; j++) {
+ C[i][j] = -A[i][j];
+ }
+ }
+ return X;
+ }
+
+ /** C = A + B
+ @param B another matrix
+ @return A + B
+ */
+
+ public Matrix plus (Matrix B) {
+ checkMatrixDimensions(B);
+ Matrix X = new Matrix(m,n);
+ double[][] C = X.getArray();
+ for (int i = 0; i < m; i++) {
+ for (int j = 0; j < n; j++) {
+ C[i][j] = A[i][j] + B.A[i][j];
+ }
+ }
+ return X;
+ }
+
+ /** A = A + B
+ @param B another matrix
+ @return A + B
+ */
+
+ public Matrix plusEquals (Matrix B) {
+ checkMatrixDimensions(B);
+ for (int i = 0; i < m; i++) {
+ for (int j = 0; j < n; j++) {
+ A[i][j] = A[i][j] + B.A[i][j];
+ }
+ }
+ return this;
+ }
+
+ /** C = A - B
+ @param B another matrix
+ @return A - B
+ */
+
+ public Matrix minus (Matrix B) {
+ checkMatrixDimensions(B);
+ Matrix X = new Matrix(m,n);
+ double[][] C = X.getArray();
+ for (int i = 0; i < m; i++) {
+ for (int j = 0; j < n; j++) {
+ C[i][j] = A[i][j] - B.A[i][j];
+ }
+ }
+ return X;
+ }
+
+ /** A = A - B
+ @param B another matrix
+ @return A - B
+ */
+
+ public Matrix minusEquals (Matrix B) {
+ checkMatrixDimensions(B);
+ for (int i = 0; i < m; i++) {
+ for (int j = 0; j < n; j++) {
+ A[i][j] = A[i][j] - B.A[i][j];
+ }
+ }
+ return this;
+ }
+
+ /** Element-by-element multiplication, C = A.*B
+ @param B another matrix
+ @return A.*B
+ */
+
+ public Matrix arrayTimes (Matrix B) {
+ checkMatrixDimensions(B);
+ Matrix X = new Matrix(m,n);
+ double[][] C = X.getArray();
+ for (int i = 0; i < m; i++) {
+ for (int j = 0; j < n; j++) {
+ C[i][j] = A[i][j] * B.A[i][j];
+ }
+ }
+ return X;
+ }
+
+ /** Element-by-element multiplication in place, A = A.*B
+ @param B another matrix
+ @return A.*B
+ */
+
+ public Matrix arrayTimesEquals (Matrix B) {
+ checkMatrixDimensions(B);
+ for (int i = 0; i < m; i++) {
+ for (int j = 0; j < n; j++) {
+ A[i][j] = A[i][j] * B.A[i][j];
+ }
+ }
+ return this;
+ }
+
+ /** Element-by-element right division, C = A./B
+ @param B another matrix
+ @return A./B
+ */
+
+ public Matrix arrayRightDivide (Matrix B) {
+ checkMatrixDimensions(B);
+ Matrix X = new Matrix(m,n);
+ double[][] C = X.getArray();
+ for (int i = 0; i < m; i++) {
+ for (int j = 0; j < n; j++) {
+ C[i][j] = A[i][j] / B.A[i][j];
+ }
+ }
+ return X;
+ }
+
+ /** Element-by-element right division in place, A = A./B
+ @param B another matrix
+ @return A./B
+ */
+
+ public Matrix arrayRightDivideEquals (Matrix B) {
+ checkMatrixDimensions(B);
+ for (int i = 0; i < m; i++) {
+ for (int j = 0; j < n; j++) {
+ A[i][j] = A[i][j] / B.A[i][j];
+ }
+ }
+ return this;
+ }
+
+ /** Element-by-element left division, C = A.\B
+ @param B another matrix
+ @return A.\B
+ */
+
+ public Matrix arrayLeftDivide (Matrix B) {
+ checkMatrixDimensions(B);
+ Matrix X = new Matrix(m,n);
+ double[][] C = X.getArray();
+ for (int i = 0; i < m; i++) {
+ for (int j = 0; j < n; j++) {
+ C[i][j] = B.A[i][j] / A[i][j];
+ }
+ }
+ return X;
+ }
+
+ /** Element-by-element left division in place, A = A.\B
+ @param B another matrix
+ @return A.\B
+ */
+
+ public Matrix arrayLeftDivideEquals (Matrix B) {
+ checkMatrixDimensions(B);
+ for (int i = 0; i < m; i++) {
+ for (int j = 0; j < n; j++) {
+ A[i][j] = B.A[i][j] / A[i][j];
+ }
+ }
+ return this;
+ }
+
+ /** Multiply a matrix by a scalar, C = s*A
+ @param s scalar
+ @return s*A
+ */
+
+ public Matrix times (double s) {
+ Matrix X = new Matrix(m,n);
+ double[][] C = X.getArray();
+ for (int i = 0; i < m; i++) {
+ for (int j = 0; j < n; j++) {
+ C[i][j] = s*A[i][j];
+ }
+ }
+ return X;
+ }
+
+ /** Multiply a matrix by a scalar in place, A = s*A
+ @param s scalar
+ @return replace A by s*A
+ */
+
+ public Matrix timesEquals (double s) {
+ for (int i = 0; i < m; i++) {
+ for (int j = 0; j < n; j++) {
+ A[i][j] = s*A[i][j];
+ }
+ }
+ return this;
+ }
+
+ /** Linear algebraic matrix multiplication, A * B
+ @param B another matrix
+ @return Matrix product, A * B
+ @exception IllegalArgumentException Matrix inner dimensions must agree.
+ */
+
+ public Matrix times (Matrix B) {
+ if (B.m != n) {
+ throw new IllegalArgumentException("Matrix inner dimensions must agree.");
+ }
+ Matrix X = new Matrix(m,B.n);
+ double[][] C = X.getArray();
+ double[] Bcolj = new double[n];
+ for (int j = 0; j < B.n; j++) {
+ for (int k = 0; k < n; k++) {
+ Bcolj[k] = B.A[k][j];
+ }
+ for (int i = 0; i < m; i++) {
+ double[] Arowi = A[i];
+ double s = 0;
+ for (int k = 0; k < n; k++) {
+ s += Arowi[k]*Bcolj[k];
+ }
+ C[i][j] = s;
+ }
+ }
+ return X;
+ }
+
+ /** LU Decomposition
+ @return LUDecomposition
+ @see LUDecomposition
+ */
+
+ public LUDecomposition lu () {
+ return new LUDecomposition(this);
+ }
+
+ /** QR Decomposition
+ @return QRDecomposition
+ @see QRDecomposition
+ */
+
+ public QRDecomposition qr () {
+ return new QRDecomposition(this);
+ }
+
+ /** Cholesky Decomposition
+ @return CholeskyDecomposition
+ @see CholeskyDecomposition
+ */
+
+ public CholeskyDecomposition chol () {
+ return new CholeskyDecomposition(this);
+ }
+
+ /** Singular Value Decomposition
+ @return SingularValueDecomposition
+ @see SingularValueDecomposition
+ */
+
+ public SingularValueDecomposition svd () {
+ return new SingularValueDecomposition(this);
+ }
+
+ /** Eigenvalue Decomposition
+ @return EigenvalueDecomposition
+ @see EigenvalueDecomposition
+ */
+
+ public EigenvalueDecomposition eig () {
+ return new EigenvalueDecomposition(this);
+ }
+
+ /** Solve A*X = B
+ @param B right hand side
+ @return solution if A is square, least squares solution otherwise
+ */
+
+ public Matrix solve (Matrix B) {
+ return (m == n ? (new LUDecomposition(this)).solve(B) :
+ (new QRDecomposition(this)).solve(B));
+ }
+
+ /** Solve X*A = B, which is also A'*X' = B'
+ @param B right hand side
+ @return solution if A is square, least squares solution otherwise.
+ */
+
+ public Matrix solveTranspose (Matrix B) {
+ return transpose().solve(B.transpose());
+ }
+
+ /** Matrix inverse or pseudoinverse
+ @return inverse(A) if A is square, pseudoinverse otherwise.
+ */
+
+ public Matrix inverse () {
+ return solve(identity(m,m));
+ }
+
+ /** Matrix determinant
+ @return determinant
+ */
+
+ public double det () {
+ return new LUDecomposition(this).det();
+ }
+
+ /** Matrix rank
+ @return effective numerical rank, obtained from SVD.
+ */
+
+ public int rank () {
+ return new SingularValueDecomposition(this).rank();
+ }
+
+ /** Matrix condition (2 norm)
+ @return ratio of largest to smallest singular value.
+ */
+
+ public double cond () {
+ return new SingularValueDecomposition(this).cond();
+ }
+
+ /** Matrix trace.
+ @return sum of the diagonal elements.
+ */
+
+ public double trace () {
+ double t = 0;
+ for (int i = 0; i < Math.min(m,n); i++) {
+ t += A[i][i];
+ }
+ return t;
+ }
+
+ /** Generate matrix with random elements
+ @param m Number of rows.
+ @param n Number of colums.
+ @return An m-by-n matrix with uniformly distributed random elements.
+ */
+
+ public static Matrix random (int m, int n) {
+ Matrix A = new Matrix(m,n);
+ double[][] X = A.getArray();
+ for (int i = 0; i < m; i++) {
+ for (int j = 0; j < n; j++) {
+ X[i][j] = Math.random();
+ }
+ }
+ return A;
+ }
+
+ /** Generate identity matrix
+ @param m Number of rows.
+ @param n Number of colums.
+ @return An m-by-n matrix with ones on the diagonal and zeros elsewhere.
+ */
+
+ public static Matrix identity (int m, int n) {
+ Matrix A = new Matrix(m,n);
+ double[][] X = A.getArray();
+ for (int i = 0; i < m; i++) {
+ for (int j = 0; j < n; j++) {
+ X[i][j] = (i == j ? 1.0 : 0.0);
+ }
+ }
+ return A;
+ }
+
+
+ /** Print the matrix to stdout. Line the elements up in columns
+ * with a Fortran-like 'Fw.d' style format.
+ @param w Column width.
+ @param d Number of digits after the decimal.
+ */
+
+ public void print (int w, int d) {
+ print(new PrintWriter(System.out,true),w,d); }
+
+ /** Print the matrix to the output stream. Line the elements up in
+ * columns with a Fortran-like 'Fw.d' style format.
+ @param output Output stream.
+ @param w Column width.
+ @param d Number of digits after the decimal.
+ */
+
+ public void print (PrintWriter output, int w, int d) {
+ DecimalFormat format = new DecimalFormat();
+ format.setDecimalFormatSymbols(new DecimalFormatSymbols(Locale.US));
+ format.setMinimumIntegerDigits(1);
+ format.setMaximumFractionDigits(d);
+ format.setMinimumFractionDigits(d);
+ format.setGroupingUsed(false);
+ print(output,format,w+2);
+ }
+
+ /** Print the matrix to stdout. Line the elements up in columns.
+ * Use the format object, and right justify within columns of width
+ * characters.
+ * Note that is the matrix is to be read back in, you probably will want
+ * to use a NumberFormat that is set to US Locale.
+ @param format A Formatting object for individual elements.
+ @param width Field width for each column.
+ @see java.text.DecimalFormat#setDecimalFormatSymbols
+ */
+
+ public void print (NumberFormat format, int width) {
+ print(new PrintWriter(System.out,true),format,width); }
+
+ // DecimalFormat is a little disappointing coming from Fortran or C's printf.
+ // Since it doesn't pad on the left, the elements will come out different
+ // widths. Consequently, we'll pass the desired column width in as an
+ // argument and do the extra padding ourselves.
+
+ /** Print the matrix to the output stream. Line the elements up in columns.
+ * Use the format object, and right justify within columns of width
+ * characters.
+ * Note that is the matrix is to be read back in, you probably will want
+ * to use a NumberFormat that is set to US Locale.
+ @param output the output stream.
+ @param format A formatting object to format the matrix elements
+ @param width Column width.
+ @see java.text.DecimalFormat#setDecimalFormatSymbols
+ */
+
+ public void print (PrintWriter output, NumberFormat format, int width) {
+ output.println(); // start on new line.
+ for (int i = 0; i < m; i++) {
+ for (int j = 0; j < n; j++) {
+ String s = format.format(A[i][j]); // format the number
+ int padding = Math.max(1,width-s.length()); // At _least_ 1 space
+ for (int k = 0; k < padding; k++)
+ output.print(' ');
+ output.print(s);
+ }
+ output.println();
+ }
+ output.println(); // end with blank line.
+ }
+
+ /** Read a matrix from a stream. The format is the same the print method,
+ * so printed matrices can be read back in (provided they were printed using
+ * US Locale). Elements are separated by
+ * whitespace, all the elements for each row appear on a single line,
+ * the last row is followed by a blank line.
+ @param input the input stream.
+ */
+
+ public static Matrix read (BufferedReader input) throws java.io.IOException {
+ StreamTokenizer tokenizer= new StreamTokenizer(input);
+
+ // Although StreamTokenizer will parse numbers, it doesn't recognize
+ // scientific notation (E or D); however, Double.valueOf does.
+ // The strategy here is to disable StreamTokenizer's number parsing.
+ // We'll only get whitespace delimited words, EOL's and EOF's.
+ // These words should all be numbers, for Double.valueOf to parse.
+
+ tokenizer.resetSyntax();
+ tokenizer.wordChars(0,255);
+ tokenizer.whitespaceChars(0, ' ');
+ tokenizer.eolIsSignificant(true);
+ java.util.Vector vD = new java.util.Vector();
+
+ // Ignore initial empty lines
+ while (tokenizer.nextToken() == StreamTokenizer.TT_EOL);
+ if (tokenizer.ttype == StreamTokenizer.TT_EOF)
+ throw new java.io.IOException("Unexpected EOF on matrix read.");
+ do {
+ vD.addElement(Double.valueOf(tokenizer.sval)); // Read & store 1st row.
+ } while (tokenizer.nextToken() == StreamTokenizer.TT_WORD);
+
+ int n = vD.size(); // Now we've got the number of columns!
+ double row[] = new double[n];
+ for (int j=0; j v = new java.util.Vector();
+ v.addElement(row); // Start storing rows instead of columns.
+ while (tokenizer.nextToken() == StreamTokenizer.TT_WORD) {
+ // While non-empty lines
+ v.addElement(row = new double[n]);
+ int j = 0;
+ do {
+ if (j >= n) throw new java.io.IOException
+ ("Row " + v.size() + " is too long.");
+ row[j++] = Double.valueOf(tokenizer.sval).doubleValue();
+ } while (tokenizer.nextToken() == StreamTokenizer.TT_WORD);
+ if (j < n) throw new java.io.IOException
+ ("Row " + v.size() + " is too short.");
+ }
+ int m = v.size(); // Now we've got the number of rows.
+ double[][] A = new double[m][];
+ v.copyInto(A); // copy the rows out of the vector
+ return new Matrix(A);
+ }
+
+
+/* ------------------------
+ Private Methods
+ * ------------------------ */
+
+ /** Check if size(A) == size(B) **/
+
+ private void checkMatrixDimensions (Matrix B) {
+ if (B.m != m || B.n != n) {
+ throw new IllegalArgumentException("Matrix dimensions must agree.");
+ }
+ }
+
+ private static final long serialVersionUID = 1;
+}
diff --git a/Java/src/Jama/QRDecomposition.java b/Java/src/Jama/QRDecomposition.java
new file mode 100644
index 00000000..10916be4
--- /dev/null
+++ b/Java/src/Jama/QRDecomposition.java
@@ -0,0 +1,219 @@
+package Jama;
+import Jama.util.*;
+
+/** QR Decomposition.
+
+ For an m-by-n matrix A with m >= n, the QR decomposition is an m-by-n
+ orthogonal matrix Q and an n-by-n upper triangular matrix R so that
+ A = Q*R.
+
+ The QR decompostion always exists, even if the matrix does not have
+ full rank, so the constructor will never fail. The primary use of the
+ QR decomposition is in the least squares solution of nonsquare systems
+ of simultaneous linear equations. This will fail if isFullRank()
+ returns false.
+*/
+
+public class QRDecomposition implements java.io.Serializable {
+
+/* ------------------------
+ Class variables
+ * ------------------------ */
+
+ /** Array for internal storage of decomposition.
+ @serial internal array storage.
+ */
+ private double[][] QR;
+
+ /** Row and column dimensions.
+ @serial column dimension.
+ @serial row dimension.
+ */
+ private int m, n;
+
+ /** Array for internal storage of diagonal of R.
+ @serial diagonal of R.
+ */
+ private double[] Rdiag;
+
+/* ------------------------
+ Constructor
+ * ------------------------ */
+
+ /** QR Decomposition, computed by Householder reflections.
+ Structure to access R and the Householder vectors and compute Q.
+ @param A Rectangular matrix
+ */
+
+ public QRDecomposition (Matrix A) {
+ // Initialize.
+ QR = A.getArrayCopy();
+ m = A.getRowDimension();
+ n = A.getColumnDimension();
+ Rdiag = new double[n];
+
+ // Main loop.
+ for (int k = 0; k < n; k++) {
+ // Compute 2-norm of k-th column without under/overflow.
+ double nrm = 0;
+ for (int i = k; i < m; i++) {
+ nrm = Maths.hypot(nrm,QR[i][k]);
+ }
+
+ if (nrm != 0.0) {
+ // Form k-th Householder vector.
+ if (QR[k][k] < 0) {
+ nrm = -nrm;
+ }
+ for (int i = k; i < m; i++) {
+ QR[i][k] /= nrm;
+ }
+ QR[k][k] += 1.0;
+
+ // Apply transformation to remaining columns.
+ for (int j = k+1; j < n; j++) {
+ double s = 0.0;
+ for (int i = k; i < m; i++) {
+ s += QR[i][k]*QR[i][j];
+ }
+ s = -s/QR[k][k];
+ for (int i = k; i < m; i++) {
+ QR[i][j] += s*QR[i][k];
+ }
+ }
+ }
+ Rdiag[k] = -nrm;
+ }
+ }
+
+/* ------------------------
+ Public Methods
+ * ------------------------ */
+
+ /** Is the matrix full rank?
+ @return true if R, and hence A, has full rank.
+ */
+
+ public boolean isFullRank () {
+ for (int j = 0; j < n; j++) {
+ if (Rdiag[j] == 0)
+ return false;
+ }
+ return true;
+ }
+
+ /** Return the Householder vectors
+ @return Lower trapezoidal matrix whose columns define the reflections
+ */
+
+ public Matrix getH () {
+ Matrix X = new Matrix(m,n);
+ double[][] H = X.getArray();
+ for (int i = 0; i < m; i++) {
+ for (int j = 0; j < n; j++) {
+ if (i >= j) {
+ H[i][j] = QR[i][j];
+ } else {
+ H[i][j] = 0.0;
+ }
+ }
+ }
+ return X;
+ }
+
+ /** Return the upper triangular factor
+ @return R
+ */
+
+ public Matrix getR () {
+ Matrix X = new Matrix(n,n);
+ double[][] R = X.getArray();
+ for (int i = 0; i < n; i++) {
+ for (int j = 0; j < n; j++) {
+ if (i < j) {
+ R[i][j] = QR[i][j];
+ } else if (i == j) {
+ R[i][j] = Rdiag[i];
+ } else {
+ R[i][j] = 0.0;
+ }
+ }
+ }
+ return X;
+ }
+
+ /** Generate and return the (economy-sized) orthogonal factor
+ @return Q
+ */
+
+ public Matrix getQ () {
+ Matrix X = new Matrix(m,n);
+ double[][] Q = X.getArray();
+ for (int k = n-1; k >= 0; k--) {
+ for (int i = 0; i < m; i++) {
+ Q[i][k] = 0.0;
+ }
+ Q[k][k] = 1.0;
+ for (int j = k; j < n; j++) {
+ if (QR[k][k] != 0) {
+ double s = 0.0;
+ for (int i = k; i < m; i++) {
+ s += QR[i][k]*Q[i][j];
+ }
+ s = -s/QR[k][k];
+ for (int i = k; i < m; i++) {
+ Q[i][j] += s*QR[i][k];
+ }
+ }
+ }
+ }
+ return X;
+ }
+
+ /** Least squares solution of A*X = B
+ @param B A Matrix with as many rows as A and any number of columns.
+ @return X that minimizes the two norm of Q*R*X-B.
+ @exception IllegalArgumentException Matrix row dimensions must agree.
+ @exception RuntimeException Matrix is rank deficient.
+ */
+
+ public Matrix solve (Matrix B) {
+ if (B.getRowDimension() != m) {
+ throw new IllegalArgumentException("Matrix row dimensions must agree.");
+ }
+ if (!this.isFullRank()) {
+ throw new RuntimeException("Matrix is rank deficient.");
+ }
+
+ // Copy right hand side
+ int nx = B.getColumnDimension();
+ double[][] X = B.getArrayCopy();
+
+ // Compute Y = transpose(Q)*B
+ for (int k = 0; k < n; k++) {
+ for (int j = 0; j < nx; j++) {
+ double s = 0.0;
+ for (int i = k; i < m; i++) {
+ s += QR[i][k]*X[i][j];
+ }
+ s = -s/QR[k][k];
+ for (int i = k; i < m; i++) {
+ X[i][j] += s*QR[i][k];
+ }
+ }
+ }
+ // Solve R*X = Y;
+ for (int k = n-1; k >= 0; k--) {
+ for (int j = 0; j < nx; j++) {
+ X[k][j] /= Rdiag[k];
+ }
+ for (int i = 0; i < k; i++) {
+ for (int j = 0; j < nx; j++) {
+ X[i][j] -= X[k][j]*QR[i][k];
+ }
+ }
+ }
+ return (new Matrix(X,n,nx).getMatrix(0,n-1,0,nx-1));
+ }
+ private static final long serialVersionUID = 1;
+}
diff --git a/Java/src/Jama/SingularValueDecomposition.java b/Java/src/Jama/SingularValueDecomposition.java
new file mode 100644
index 00000000..b289e82f
--- /dev/null
+++ b/Java/src/Jama/SingularValueDecomposition.java
@@ -0,0 +1,548 @@
+package Jama;
+import Jama.util.*;
+
+ /** Singular Value Decomposition.
+
+ For an m-by-n matrix A with m >= n, the singular value decomposition is
+ an m-by-n orthogonal matrix U, an n-by-n diagonal matrix S, and
+ an n-by-n orthogonal matrix V so that A = U*S*V'.
+
+ The singular values, sigma[k] = S[k][k], are ordered so that
+ sigma[0] >= sigma[1] >= ... >= sigma[n-1].
+
+ The singular value decompostion always exists, so the constructor will
+ never fail. The matrix condition number and the effective numerical
+ rank can be computed from this decomposition.
+ */
+
+public class SingularValueDecomposition implements java.io.Serializable {
+
+/* ------------------------
+ Class variables
+ * ------------------------ */
+
+ /** Arrays for internal storage of U and V.
+ @serial internal storage of U.
+ @serial internal storage of V.
+ */
+ private double[][] U, V;
+
+ /** Array for internal storage of singular values.
+ @serial internal storage of singular values.
+ */
+ private double[] s;
+
+ /** Row and column dimensions.
+ @serial row dimension.
+ @serial column dimension.
+ */
+ private int m, n;
+
+/* ------------------------
+ Constructor
+ * ------------------------ */
+
+ /** Construct the singular value decomposition
+ Structure to access U, S and V.
+ @param Arg Rectangular matrix
+ */
+
+ public SingularValueDecomposition (Matrix Arg) {
+
+ // Derived from LINPACK code.
+ // Initialize.
+ double[][] A = Arg.getArrayCopy();
+ m = Arg.getRowDimension();
+ n = Arg.getColumnDimension();
+
+ /* Apparently the failing cases are only a proper subset of (m= n"); }
+ */
+ int nu = Math.min(m,n);
+ s = new double [Math.min(m+1,n)];
+ U = new double [m][nu];
+ V = new double [n][n];
+ double[] e = new double [n];
+ double[] work = new double [m];
+ boolean wantu = true;
+ boolean wantv = true;
+
+ // Reduce A to bidiagonal form, storing the diagonal elements
+ // in s and the super-diagonal elements in e.
+
+ int nct = Math.min(m-1,n);
+ int nrt = Math.max(0,Math.min(n-2,m));
+ for (int k = 0; k < Math.max(nct,nrt); k++) {
+ if (k < nct) {
+
+ // Compute the transformation for the k-th column and
+ // place the k-th diagonal in s[k].
+ // Compute 2-norm of k-th column without under/overflow.
+ s[k] = 0;
+ for (int i = k; i < m; i++) {
+ s[k] = Maths.hypot(s[k],A[i][k]);
+ }
+ if (s[k] != 0.0) {
+ if (A[k][k] < 0.0) {
+ s[k] = -s[k];
+ }
+ for (int i = k; i < m; i++) {
+ A[i][k] /= s[k];
+ }
+ A[k][k] += 1.0;
+ }
+ s[k] = -s[k];
+ }
+ for (int j = k+1; j < n; j++) {
+ if ((k < nct) & (s[k] != 0.0)) {
+
+ // Apply the transformation.
+
+ double t = 0;
+ for (int i = k; i < m; i++) {
+ t += A[i][k]*A[i][j];
+ }
+ t = -t/A[k][k];
+ for (int i = k; i < m; i++) {
+ A[i][j] += t*A[i][k];
+ }
+ }
+
+ // Place the k-th row of A into e for the
+ // subsequent calculation of the row transformation.
+
+ e[j] = A[k][j];
+ }
+ if (wantu & (k < nct)) {
+
+ // Place the transformation in U for subsequent back
+ // multiplication.
+
+ for (int i = k; i < m; i++) {
+ U[i][k] = A[i][k];
+ }
+ }
+ if (k < nrt) {
+
+ // Compute the k-th row transformation and place the
+ // k-th super-diagonal in e[k].
+ // Compute 2-norm without under/overflow.
+ e[k] = 0;
+ for (int i = k+1; i < n; i++) {
+ e[k] = Maths.hypot(e[k],e[i]);
+ }
+ if (e[k] != 0.0) {
+ if (e[k+1] < 0.0) {
+ e[k] = -e[k];
+ }
+ for (int i = k+1; i < n; i++) {
+ e[i] /= e[k];
+ }
+ e[k+1] += 1.0;
+ }
+ e[k] = -e[k];
+ if ((k+1 < m) & (e[k] != 0.0)) {
+
+ // Apply the transformation.
+
+ for (int i = k+1; i < m; i++) {
+ work[i] = 0.0;
+ }
+ for (int j = k+1; j < n; j++) {
+ for (int i = k+1; i < m; i++) {
+ work[i] += e[j]*A[i][j];
+ }
+ }
+ for (int j = k+1; j < n; j++) {
+ double t = -e[j]/e[k+1];
+ for (int i = k+1; i < m; i++) {
+ A[i][j] += t*work[i];
+ }
+ }
+ }
+ if (wantv) {
+
+ // Place the transformation in V for subsequent
+ // back multiplication.
+
+ for (int i = k+1; i < n; i++) {
+ V[i][k] = e[i];
+ }
+ }
+ }
+ }
+
+ // Set up the final bidiagonal matrix or order p.
+
+ int p = Math.min(n,m+1);
+ if (nct < n) {
+ s[nct] = A[nct][nct];
+ }
+ if (m < p) {
+ s[p-1] = 0.0;
+ }
+ if (nrt+1 < p) {
+ e[nrt] = A[nrt][p-1];
+ }
+ e[p-1] = 0.0;
+
+ // If required, generate U.
+
+ if (wantu) {
+ for (int j = nct; j < nu; j++) {
+ for (int i = 0; i < m; i++) {
+ U[i][j] = 0.0;
+ }
+ U[j][j] = 1.0;
+ }
+ for (int k = nct-1; k >= 0; k--) {
+ if (s[k] != 0.0) {
+ for (int j = k+1; j < nu; j++) {
+ double t = 0;
+ for (int i = k; i < m; i++) {
+ t += U[i][k]*U[i][j];
+ }
+ t = -t/U[k][k];
+ for (int i = k; i < m; i++) {
+ U[i][j] += t*U[i][k];
+ }
+ }
+ for (int i = k; i < m; i++ ) {
+ U[i][k] = -U[i][k];
+ }
+ U[k][k] = 1.0 + U[k][k];
+ for (int i = 0; i < k-1; i++) {
+ U[i][k] = 0.0;
+ }
+ } else {
+ for (int i = 0; i < m; i++) {
+ U[i][k] = 0.0;
+ }
+ U[k][k] = 1.0;
+ }
+ }
+ }
+
+ // If required, generate V.
+
+ if (wantv) {
+ for (int k = n-1; k >= 0; k--) {
+ if ((k < nrt) & (e[k] != 0.0)) {
+ for (int j = k+1; j < nu; j++) {
+ double t = 0;
+ for (int i = k+1; i < n; i++) {
+ t += V[i][k]*V[i][j];
+ }
+ t = -t/V[k+1][k];
+ for (int i = k+1; i < n; i++) {
+ V[i][j] += t*V[i][k];
+ }
+ }
+ }
+ for (int i = 0; i < n; i++) {
+ V[i][k] = 0.0;
+ }
+ V[k][k] = 1.0;
+ }
+ }
+
+ // Main iteration loop for the singular values.
+
+ int pp = p-1;
+ int iter = 0;
+ double eps = Math.pow(2.0,-52.0);
+ double tiny = Math.pow(2.0,-966.0);
+ while (p > 0) {
+ int k,kase;
+
+ // Here is where a test for too many iterations would go.
+
+ // This section of the program inspects for
+ // negligible elements in the s and e arrays. On
+ // completion the variables kase and k are set as follows.
+
+ // kase = 1 if s(p) and e[k-1] are negligible and k
= -1; k--) {
+ if (k == -1) {
+ break;
+ }
+ if (Math.abs(e[k]) <=
+ tiny + eps*(Math.abs(s[k]) + Math.abs(s[k+1]))) {
+ e[k] = 0.0;
+ break;
+ }
+ }
+ if (k == p-2) {
+ kase = 4;
+ } else {
+ int ks;
+ for (ks = p-1; ks >= k; ks--) {
+ if (ks == k) {
+ break;
+ }
+ double t = (ks != p ? Math.abs(e[ks]) : 0.) +
+ (ks != k+1 ? Math.abs(e[ks-1]) : 0.);
+ if (Math.abs(s[ks]) <= tiny + eps*t) {
+ s[ks] = 0.0;
+ break;
+ }
+ }
+ if (ks == k) {
+ kase = 3;
+ } else if (ks == p-1) {
+ kase = 1;
+ } else {
+ kase = 2;
+ k = ks;
+ }
+ }
+ k++;
+
+ // Perform the task indicated by kase.
+
+ switch (kase) {
+
+ // Deflate negligible s(p).
+
+ case 1: {
+ double f = e[p-2];
+ e[p-2] = 0.0;
+ for (int j = p-2; j >= k; j--) {
+ double t = Maths.hypot(s[j],f);
+ double cs = s[j]/t;
+ double sn = f/t;
+ s[j] = t;
+ if (j != k) {
+ f = -sn*e[j-1];
+ e[j-1] = cs*e[j-1];
+ }
+ if (wantv) {
+ for (int i = 0; i < n; i++) {
+ t = cs*V[i][j] + sn*V[i][p-1];
+ V[i][p-1] = -sn*V[i][j] + cs*V[i][p-1];
+ V[i][j] = t;
+ }
+ }
+ }
+ }
+ break;
+
+ // Split at negligible s(k).
+
+ case 2: {
+ double f = e[k-1];
+ e[k-1] = 0.0;
+ for (int j = k; j < p; j++) {
+ double t = Maths.hypot(s[j],f);
+ double cs = s[j]/t;
+ double sn = f/t;
+ s[j] = t;
+ f = -sn*e[j];
+ e[j] = cs*e[j];
+ if (wantu) {
+ for (int i = 0; i < m; i++) {
+ t = cs*U[i][j] + sn*U[i][k-1];
+ U[i][k-1] = -sn*U[i][j] + cs*U[i][k-1];
+ U[i][j] = t;
+ }
+ }
+ }
+ }
+ break;
+
+ // Perform one qr step.
+
+ case 3: {
+
+ // Calculate the shift.
+
+ double scale = Math.max(Math.max(Math.max(Math.max(
+ Math.abs(s[p-1]),Math.abs(s[p-2])),Math.abs(e[p-2])),
+ Math.abs(s[k])),Math.abs(e[k]));
+ double sp = s[p-1]/scale;
+ double spm1 = s[p-2]/scale;
+ double epm1 = e[p-2]/scale;
+ double sk = s[k]/scale;
+ double ek = e[k]/scale;
+ double b = ((spm1 + sp)*(spm1 - sp) + epm1*epm1)/2.0;
+ double c = (sp*epm1)*(sp*epm1);
+ double shift = 0.0;
+ if ((b != 0.0) | (c != 0.0)) {
+ shift = Math.sqrt(b*b + c);
+ if (b < 0.0) {
+ shift = -shift;
+ }
+ shift = c/(b + shift);
+ }
+ double f = (sk + sp)*(sk - sp) + shift;
+ double g = sk*ek;
+
+ // Chase zeros.
+
+ for (int j = k; j < p-1; j++) {
+ double t = Maths.hypot(f,g);
+ double cs = f/t;
+ double sn = g/t;
+ if (j != k) {
+ e[j-1] = t;
+ }
+ f = cs*s[j] + sn*e[j];
+ e[j] = cs*e[j] - sn*s[j];
+ g = sn*s[j+1];
+ s[j+1] = cs*s[j+1];
+ if (wantv) {
+ for (int i = 0; i < n; i++) {
+ t = cs*V[i][j] + sn*V[i][j+1];
+ V[i][j+1] = -sn*V[i][j] + cs*V[i][j+1];
+ V[i][j] = t;
+ }
+ }
+ t = Maths.hypot(f,g);
+ cs = f/t;
+ sn = g/t;
+ s[j] = t;
+ f = cs*e[j] + sn*s[j+1];
+ s[j+1] = -sn*e[j] + cs*s[j+1];
+ g = sn*e[j+1];
+ e[j+1] = cs*e[j+1];
+ if (wantu && (j < m-1)) {
+ for (int i = 0; i < m; i++) {
+ t = cs*U[i][j] + sn*U[i][j+1];
+ U[i][j+1] = -sn*U[i][j] + cs*U[i][j+1];
+ U[i][j] = t;
+ }
+ }
+ }
+ e[p-2] = f;
+ iter = iter + 1;
+ }
+ break;
+
+ // Convergence.
+
+ case 4: {
+
+ // Make the singular values positive.
+
+ if (s[k] <= 0.0) {
+ s[k] = (s[k] < 0.0 ? -s[k] : 0.0);
+ if (wantv) {
+ for (int i = 0; i <= pp; i++) {
+ V[i][k] = -V[i][k];
+ }
+ }
+ }
+
+ // Order the singular values.
+
+ while (k < pp) {
+ if (s[k] >= s[k+1]) {
+ break;
+ }
+ double t = s[k];
+ s[k] = s[k+1];
+ s[k+1] = t;
+ if (wantv && (k < n-1)) {
+ for (int i = 0; i < n; i++) {
+ t = V[i][k+1]; V[i][k+1] = V[i][k]; V[i][k] = t;
+ }
+ }
+ if (wantu && (k < m-1)) {
+ for (int i = 0; i < m; i++) {
+ t = U[i][k+1]; U[i][k+1] = U[i][k]; U[i][k] = t;
+ }
+ }
+ k++;
+ }
+ iter = 0;
+ p--;
+ }
+ break;
+ }
+ }
+ }
+
+/* ------------------------
+ Public Methods
+ * ------------------------ */
+
+ /** Return the left singular vectors
+ @return U
+ */
+
+ public Matrix getU () {
+ return new Matrix(U,m,Math.min(m+1,n));
+ }
+
+ /** Return the right singular vectors
+ @return V
+ */
+
+ public Matrix getV () {
+ return new Matrix(V,n,n);
+ }
+
+ /** Return the one-dimensional array of singular values
+ @return diagonal of S.
+ */
+
+ public double[] getSingularValues () {
+ return s;
+ }
+
+ /** Return the diagonal matrix of singular values
+ @return S
+ */
+
+ public Matrix getS () {
+ Matrix X = new Matrix(n,n);
+ double[][] S = X.getArray();
+ for (int i = 0; i < n; i++) {
+ for (int j = 0; j < n; j++) {
+ S[i][j] = 0.0;
+ }
+ S[i][i] = this.s[i];
+ }
+ return X;
+ }
+
+ /** Two norm
+ @return max(S)
+ */
+
+ public double norm2 () {
+ return s[0];
+ }
+
+ /** Two norm condition number
+ @return max(S)/min(S)
+ */
+
+ public double cond () {
+ return s[0]/s[Math.min(m,n)-1];
+ }
+
+ /** Effective numerical matrix rank
+ @return Number of nonnegligible singular values.
+ */
+
+ public int rank () {
+ double eps = Math.pow(2.0,-52.0);
+ double tol = Math.max(m,n)*s[0]*eps;
+ int r = 0;
+ for (int i = 0; i < s.length; i++) {
+ if (s[i] > tol) {
+ r++;
+ }
+ }
+ return r;
+ }
+ private static final long serialVersionUID = 1;
+}
diff --git a/Java/src/Jama/doc/Jama/CholeskyDecomposition.html b/Java/src/Jama/doc/Jama/CholeskyDecomposition.html
new file mode 100644
index 00000000..4f1721c2
--- /dev/null
+++ b/Java/src/Jama/doc/Jama/CholeskyDecomposition.html
@@ -0,0 +1,328 @@
+
+
+
+
+ For a symmetric, positive definite matrix A, the Cholesky decomposition
+ is an lower triangular matrix L so that A = L*L'.
+
+ If the matrix is not symmetric or positive definite, the constructor
+ returns a partial decomposition and sets an internal flag that may
+ be queried by the isSPD() method.
+
+ If A is symmetric, then A = V*D*V' where the eigenvalue matrix D is
+ diagonal and the eigenvector matrix V is orthogonal.
+ I.e. A = V.times(D.times(V.transpose())) and
+ V.times(V.transpose()) equals the identity matrix.
+
+ If A is not symmetric, then the eigenvalue matrix D is block diagonal
+ with the real eigenvalues in 1-by-1 blocks and any complex eigenvalues,
+ lambda + i*mu, in 2-by-2 blocks, [lambda, mu; -mu, lambda]. The
+ columns of V represent the eigenvectors in the sense that A*V = V*D,
+ i.e. A.times(V) equals V.times(D). The matrix V may be badly
+ conditioned, or even singular, so the validity of the equation
+ A = V*D*inverse(V) depends upon V.cond().
+
+ For an m-by-n matrix A with m >= n, the LU decomposition is an m-by-n
+ unit lower triangular matrix L, an n-by-n upper triangular matrix U,
+ and a permutation vector piv of length m so that A(piv,:) = L*U.
+ If m < n, then L is m-by-m and U is m-by-n.
+
+ The LU decompostion with pivoting always exists, even if the matrix is
+ singular, so the constructor will never fail. The primary use of the
+ LU decomposition is in the solution of square systems of simultaneous
+ linear equations. This will fail if isNonsingular() returns false.
+
+ The Java Matrix Class provides the fundamental operations of numerical
+ linear algebra. Various constructors create Matrices from two dimensional
+ arrays of double precision floating point numbers. Various "gets" and
+ "sets" provide access to submatrices and matrix elements. Several methods
+ implement basic matrix arithmetic, including matrix addition and
+ multiplication, matrix norms, and element-by-element array operations.
+ Methods for reading and printing matrices are also included. All the
+ operations in this version of the Matrix Class involve real matrices.
+ Complex matrices may be handled in a future version.
+
+ Five fundamental matrix decompositions, which consist of pairs or triples
+ of matrices, permutation vectors, and the like, produce results in five
+ decomposition classes. These decompositions are accessed by the Matrix
+ class to compute solutions of simultaneous linear equations, determinants,
+ inverses and other matrix functions. The five decompositions are:
+
+
Cholesky Decomposition of symmetric, positive definite matrices.
+
LU Decomposition of rectangular matrices.
+
QR Decomposition of rectangular matrices.
+
Singular Value Decomposition of rectangular matrices.
+
Eigenvalue Decomposition of both symmetric and nonsymmetric square matrices.
+
+
+
Example of use:
+
+
Solve a linear system A x = b and compute the residual norm, ||b - A x||.
+
+ double[][] vals = {{1.,2.,3},{4.,5.,6.},{7.,8.,10.}};
+ Matrix A = new Matrix(vals);
+ Matrix b = Matrix.random(3,1);
+ Matrix x = A.solve(b);
+ Matrix r = A.times(x).minus(b);
+ double rnorm = r.normInf();
+
An m-by-n matrix with ones on the diagonal and zeros elsewhere.
+
+
+
+
+
+print
+
+public void print(int w,
+ int d)
+
+
Print the matrix to stdout. Line the elements up in columns
+ with a Fortran-like 'Fw.d' style format.
+
+
+
+
+
+
Parameters:
w - Column width.
d - Number of digits after the decimal.
+
+
+
+
+
+print
+
+public void print(PrintWriter output,
+ int w,
+ int d)
+
+
Print the matrix to the output stream. Line the elements up in
+ columns with a Fortran-like 'Fw.d' style format.
+
+
+
+
+
+
Parameters:
output - Output stream.
w - Column width.
d - Number of digits after the decimal.
+
+
+
+
+
+print
+
+public void print(NumberFormat format,
+ int width)
+
+
Print the matrix to stdout. Line the elements up in columns.
+ Use the format object, and right justify within columns of width
+ characters.
+ Note that is the matrix is to be read back in, you probably will want
+ to use a NumberFormat that is set to US Locale.
+
+
+
+
+
+
Parameters:
format - A Formatting object for individual elements.
Print the matrix to the output stream. Line the elements up in columns.
+ Use the format object, and right justify within columns of width
+ characters.
+ Note that is the matrix is to be read back in, you probably will want
+ to use a NumberFormat that is set to US Locale.
+
+
+
+
+
+
Parameters:
output - the output stream.
format - A formatting object to format the matrix elements
Read a matrix from a stream. The format is the same the print method,
+ so printed matrices can be read back in (provided they were printed using
+ US Locale). Elements are separated by
+ whitespace, all the elements for each row appear on a single line,
+ the last row is followed by a blank line.
+
+ For an m-by-n matrix A with m >= n, the QR decomposition is an m-by-n
+ orthogonal matrix Q and an n-by-n upper triangular matrix R so that
+ A = Q*R.
+
+ The QR decompostion always exists, even if the matrix does not have
+ full rank, so the constructor will never fail. The primary use of the
+ QR decomposition is in the least squares solution of nonsquare systems
+ of simultaneous linear equations. This will fail if isFullRank()
+ returns false.
+
+ For an m-by-n matrix A with m >= n, the singular value decomposition is
+ an m-by-n orthogonal matrix U, an n-by-n diagonal matrix S, and
+ an n-by-n orthogonal matrix V so that A = U*S*V'.
+
+ The singular values, sigma[k] = S[k][k], are ordered so that
+ sigma[0] >= sigma[1] >= ... >= sigma[n-1].
+
+ The singular value decompostion always exists, so the constructor will
+ never fail. The matrix condition number and the effective numerical
+ rank can be computed from this decomposition.
+
+This API (Application Programming Interface) document has pages corresponding to the items in the navigation bar, described as follows.
+Package
+
+
+
+Each package has a page that contains a list of its classes and interfaces, with a summary for each. This page can contain four categories:
+
Interfaces (italic)
Classes
Enums
Exceptions
Errors
Annotation Types
+
+
+Class/Interface
+
+
+
+Each class, interface, nested class and nested interface has its own separate page. Each of these pages has three sections consisting of a class/interface description, summary tables, and detailed member descriptions:
+
Class inheritance diagram
Direct Subclasses
All Known Subinterfaces
All Known Implementing Classes
Class/interface declaration
Class/interface description
+
+
Nested Class Summary
Field Summary
Constructor Summary
Method Summary
+
+
Field Detail
Constructor Detail
Method Detail
+Each summary entry contains the first sentence from the detailed description for that item. The summary entries are alphabetical, while the detailed descriptions are in the order they appear in the source code. This preserves the logical groupings established by the programmer.
+
+
+Annotation Type
+
+
+
+Each annotation type has its own separate page with the following sections:
+
Annotation Type declaration
Annotation Type description
Required Element Summary
Optional Element Summary
Element Detail
+
+
+
+Enum
+
+
+
+Each enum has its own separate page with the following sections:
+
Enum declaration
Enum description
Enum Constant Summary
Enum Constant Detail
+
+
+Tree (Class Hierarchy)
+
+There is a Class Hierarchy page for all packages, plus a hierarchy for each package. Each hierarchy page contains a list of classes and a list of interfaces. The classes are organized by inheritance structure starting with java.lang.Object. The interfaces do not inherit from java.lang.Object.
+
When viewing the Overview page, clicking on "Tree" displays the hierarchy for all packages.
When viewing a particular package, class or interface page, clicking "Tree" displays the hierarchy for only that package.
+
+
+Deprecated API
+
+The Deprecated API page lists all of the API that have been deprecated. A deprecated API is not recommended for use, generally due to improvements, and a replacement API is usually given. Deprecated APIs may be removed in future implementations.
+
+Index
+
+The Index contains an alphabetic list of all classes, interfaces, constructors, methods, and fields.
+
+Prev/Next
+These links take you to the next or previous class, interface, package, or related page.
+Frames/No Frames
+These links show and hide the HTML frames. All pages are available with or without frames.
+
+
+Serialized Form
+Each serializable or externalizable class has a description of its serialization fields and methods. This information is of interest to re-implementors, not to developers using the API. While there is no link in the navigation bar, you can get to this information by going to any serialized class and clicking "Serialized Form" in the "See also" section of the class description.
+
+
+
+
+
+
+
+
diff --git a/Java/src/Jama/doc/stylesheet.css b/Java/src/Jama/doc/stylesheet.css
new file mode 100644
index 00000000..6ea9e516
--- /dev/null
+++ b/Java/src/Jama/doc/stylesheet.css
@@ -0,0 +1,29 @@
+/* Javadoc style sheet */
+
+/* Define colors, fonts and other style attributes here to override the defaults */
+
+/* Page background color */
+body { background-color: #FFFFFF; color:#000000 }
+
+/* Headings */
+h1 { font-size: 145% }
+
+/* Table colors */
+.TableHeadingColor { background: #CCCCFF; color:#000000 } /* Dark mauve */
+.TableSubHeadingColor { background: #EEEEFF; color:#000000 } /* Light mauve */
+.TableRowColor { background: #FFFFFF; color:#000000 } /* White */
+
+/* Font used in left-hand frame lists */
+.FrameTitleFont { font-size: 100%; font-family: Helvetica, Arial, sans-serif; color:#000000 }
+.FrameHeadingFont { font-size: 90%; font-family: Helvetica, Arial, sans-serif; color:#000000 }
+.FrameItemFont { font-size: 90%; font-family: Helvetica, Arial, sans-serif; color:#000000 }
+
+/* Navigation bar fonts and colors */
+.NavBarCell1 { background-color:#EEEEFF; color:#000000} /* Light mauve */
+.NavBarCell1Rev { background-color:#00008B; color:#FFFFFF} /* Dark Blue */
+.NavBarFont1 { font-family: Arial, Helvetica, sans-serif; color:#000000;color:#000000;}
+.NavBarFont1Rev { font-family: Arial, Helvetica, sans-serif; color:#FFFFFF;color:#FFFFFF;}
+
+.NavBarCell2 { font-family: Arial, Helvetica, sans-serif; background-color:#FFFFFF; color:#000000}
+.NavBarCell3 { font-family: Arial, Helvetica, sans-serif; background-color:#FFFFFF; color:#000000}
+
diff --git a/Java/src/Jama/examples/MagicSquareExample.java b/Java/src/Jama/examples/MagicSquareExample.java
new file mode 100644
index 00000000..cf0b7db9
--- /dev/null
+++ b/Java/src/Jama/examples/MagicSquareExample.java
@@ -0,0 +1,163 @@
+package Jama.examples;
+import Jama.*;
+import java.util.Date;
+
+/** Example of use of Matrix Class, featuring magic squares. **/
+
+public class MagicSquareExample {
+
+ /** Generate magic square test matrix. **/
+
+ public static Matrix magic (int n) {
+
+ double[][] M = new double[n][n];
+
+ // Odd order
+
+ if ((n % 2) == 1) {
+ int a = (n+1)/2;
+ int b = (n+1);
+ for (int j = 0; j < n; j++) {
+ for (int i = 0; i < n; i++) {
+ M[i][j] = n*((i+j+a) % n) + ((i+2*j+b) % n) + 1;
+ }
+ }
+
+ // Doubly Even Order
+
+ } else if ((n % 4) == 0) {
+ for (int j = 0; j < n; j++) {
+ for (int i = 0; i < n; i++) {
+ if (((i+1)/2)%2 == ((j+1)/2)%2) {
+ M[i][j] = n*n-n*i-j;
+ } else {
+ M[i][j] = n*i+j+1;
+ }
+ }
+ }
+
+ // Singly Even Order
+
+ } else {
+ int p = n/2;
+ int k = (n-2)/4;
+ Matrix A = magic(p);
+ for (int j = 0; j < p; j++) {
+ for (int i = 0; i < p; i++) {
+ double aij = A.get(i,j);
+ M[i][j] = aij;
+ M[i][j+p] = aij + 2*p*p;
+ M[i+p][j] = aij + 3*p*p;
+ M[i+p][j+p] = aij + p*p;
+ }
+ }
+ for (int i = 0; i < p; i++) {
+ for (int j = 0; j < k; j++) {
+ double t = M[i][j]; M[i][j] = M[i+p][j]; M[i+p][j] = t;
+ }
+ for (int j = n-k+1; j < n; j++) {
+ double t = M[i][j]; M[i][j] = M[i+p][j]; M[i+p][j] = t;
+ }
+ }
+ double t = M[k][0]; M[k][0] = M[k+p][0]; M[k+p][0] = t;
+ t = M[k][k]; M[k][k] = M[k+p][k]; M[k+p][k] = t;
+ }
+ return new Matrix(M);
+ }
+
+ /** Shorten spelling of print. **/
+
+ private static void print (String s) {
+ System.out.print(s);
+ }
+
+ /** Format double with Fw.d. **/
+
+ public static String fixedWidthDoubletoString (double x, int w, int d) {
+ java.text.DecimalFormat fmt = new java.text.DecimalFormat();
+ fmt.setMaximumFractionDigits(d);
+ fmt.setMinimumFractionDigits(d);
+ fmt.setGroupingUsed(false);
+ String s = fmt.format(x);
+ while (s.length() < w) {
+ s = " " + s;
+ }
+ return s;
+ }
+
+ /** Format integer with Iw. **/
+
+ public static String fixedWidthIntegertoString (int n, int w) {
+ String s = Integer.toString(n);
+ while (s.length() < w) {
+ s = " " + s;
+ }
+ return s;
+ }
+
+
+ public static void main (String argv[]) {
+
+ /*
+ | Tests LU, QR, SVD and symmetric Eig decompositions.
+ |
+ | n = order of magic square.
+ | trace = diagonal sum, should be the magic sum, (n^3 + n)/2.
+ | max_eig = maximum eigenvalue of (A + A')/2, should equal trace.
+ | rank = linear algebraic rank,
+ | should equal n if n is odd, be less than n if n is even.
+ | cond = L_2 condition number, ratio of singular values.
+ | lu_res = test of LU factorization, norm1(L*U-A(p,:))/(n*eps).
+ | qr_res = test of QR factorization, norm1(Q*R-A)/(n*eps).
+ */
+
+ print("\n Test of Matrix Class, using magic squares.\n");
+ print(" See MagicSquareExample.main() for an explanation.\n");
+ print("\n n trace max_eig rank cond lu_res qr_res\n\n");
+
+ Date start_time = new Date();
+ double eps = Math.pow(2.0,-52.0);
+ for (int n = 3; n <= 32; n++) {
+ print(fixedWidthIntegertoString(n,7));
+
+ Matrix M = magic(n);
+
+ int t = (int) M.trace();
+ print(fixedWidthIntegertoString(t,10));
+
+ EigenvalueDecomposition E =
+ new EigenvalueDecomposition(M.plus(M.transpose()).times(0.5));
+ double[] d = E.getRealEigenvalues();
+ print(fixedWidthDoubletoString(d[n-1],14,3));
+
+ int r = M.rank();
+ print(fixedWidthIntegertoString(r,7));
+
+ double c = M.cond();
+ print(c < 1/eps ? fixedWidthDoubletoString(c,12,3) :
+ " Inf");
+
+ LUDecomposition LU = new LUDecomposition(M);
+ Matrix L = LU.getL();
+ Matrix U = LU.getU();
+ int[] p = LU.getPivot();
+ Matrix R = L.times(U).minus(M.getMatrix(p,0,n-1));
+ double res = R.norm1()/(n*eps);
+ print(fixedWidthDoubletoString(res,12,3));
+
+ QRDecomposition QR = new QRDecomposition(M);
+ Matrix Q = QR.getQ();
+ R = QR.getR();
+ R = Q.times(R).minus(M);
+ res = R.norm1()/(n*eps);
+ print(fixedWidthDoubletoString(res,12,3));
+
+ print("\n");
+ }
+ Date stop_time = new Date();
+ double etime = (stop_time.getTime() - start_time.getTime())/1000.;
+ print("\nElapsed Time = " +
+ fixedWidthDoubletoString(etime,12,3) + " seconds\n");
+ print("Adios\n");
+ }
+}
diff --git a/Java/src/Jama/test/TestMatrix.java b/Java/src/Jama/test/TestMatrix.java
new file mode 100644
index 00000000..d08ca2c4
--- /dev/null
+++ b/Java/src/Jama/test/TestMatrix.java
@@ -0,0 +1,957 @@
+package Jama.test;
+import Jama.*;
+import java.io.*;
+import java.util.zip.GZIPInputStream;
+import java.text.DecimalFormat;
+import java.text.DecimalFormatSymbols;
+import java.util.Locale;
+
+/** TestMatrix tests the functionality of the Jama Matrix class and associated decompositions.
+
+Run the test from the command line using
+
+ java Jama.test.TestMatrix
+
+Detailed output is provided indicating the functionality being tested
+and whether the functionality is correctly implemented. Exception handling
+is also tested.
+
+The test is designed to run to completion and give a summary of any implementation errors
+encountered. The final output should be:
+
+ TestMatrix completed.
+ Total errors reported: n1
+ Total warning reported: n2
+
+If the test does not run to completion, this indicates that there is a
+substantial problem within the implementation that was not anticipated in the test design.
+The stopping point should give an indication of where the problem exists.
+**/
+public class TestMatrix {
+ public static void main (String argv[]) {
+ Matrix A,B,C,Z,O,I,R,S,X,SUB,M,T,SQ,DEF,SOL;
+ // Uncomment this to test IO in a different locale.
+ // Locale.setDefault(Locale.GERMAN);
+ int errorCount=0;
+ int warningCount=0;
+ double tmp, s;
+ double[] columnwise = {1.,2.,3.,4.,5.,6.,7.,8.,9.,10.,11.,12.};
+ double[] rowwise = {1.,4.,7.,10.,2.,5.,8.,11.,3.,6.,9.,12.};
+ double[][] avals = {{1.,4.,7.,10.},{2.,5.,8.,11.},{3.,6.,9.,12.}};
+ double[][] rankdef = avals;
+ double[][] tvals = {{1.,2.,3.},{4.,5.,6.},{7.,8.,9.},{10.,11.,12.}};
+ double[][] subavals = {{5.,8.,11.},{6.,9.,12.}};
+ double[][] rvals = {{1.,4.,7.},{2.,5.,8.,11.},{3.,6.,9.,12.}};
+ double[][] pvals = {{4.,1.,1.},{1.,2.,3.},{1.,3.,6.}};
+ double[][] ivals = {{1.,0.,0.,0.},{0.,1.,0.,0.},{0.,0.,1.,0.}};
+ double[][] evals =
+ {{0.,1.,0.,0.},{1.,0.,2.e-7,0.},{0.,-2.e-7,0.,1.},{0.,0.,1.,0.}};
+ double[][] square = {{166.,188.,210.},{188.,214.,240.},{210.,240.,270.}};
+ double[][] sqSolution = {{13.},{15.}};
+ double[][] condmat = {{1.,3.},{7.,9.}};
+ double[][] badeigs = {{0,0,0,0,0}, {0,0,0,0,1},{0,0,0,1,0},
+ {1,1,0,0,1},{1,0,1,0,1}};
+ int rows=3,cols=4;
+ int invalidld=5;/* should trigger bad shape for construction with val */
+ int raggedr=0; /* (raggedr,raggedc) should be out of bounds in ragged array */
+ int raggedc=4;
+ int validld=3; /* leading dimension of intended test Matrices */
+ int nonconformld=4; /* leading dimension which is valid, but nonconforming */
+ int ib=1,ie=2,jb=1,je=3; /* index ranges for sub Matrix */
+ int[] rowindexset = {1,2};
+ int[] badrowindexset = {1,3};
+ int[] columnindexset = {1,2,3};
+ int[] badcolumnindexset = {1,2,4};
+ double columnsummax = 33.;
+ double rowsummax = 30.;
+ double sumofdiagonals = 15;
+ double sumofsquares = 650;
+
+/**
+ Constructors and constructor-like methods:
+ double[], int
+ double[][]
+ int, int
+ int, int, double
+ int, int, double[][]
+ constructWithCopy(double[][])
+ random(int,int)
+ identity(int)
+**/
+
+ print("\nTesting constructors and constructor-like methods...\n");
+ try{
+ /** check that exception is thrown in packed constructor with invalid length **/
+ A = new Matrix(columnwise,invalidld);
+ errorCount = try_failure(errorCount,"Catch invalid length in packed constructor... ",
+ "exception not thrown for invalid input");
+ } catch ( IllegalArgumentException e ) {
+ try_success("Catch invalid length in packed constructor... ",
+ e.getMessage());
+ }
+ try{
+ /** check that exception is thrown in default constructor
+ if input array is 'ragged' **/
+ A = new Matrix(rvals);
+ tmp = A.get(raggedr,raggedc);
+ } catch ( IllegalArgumentException e ) {
+ try_success("Catch ragged input to default constructor... ",
+ e.getMessage());
+ } catch ( java.lang.ArrayIndexOutOfBoundsException e ) {
+ errorCount = try_failure(errorCount,"Catch ragged input to constructor... ",
+ "exception not thrown in construction...ArrayIndexOutOfBoundsException thrown later");
+ }
+ try{
+ /** check that exception is thrown in constructWithCopy
+ if input array is 'ragged' **/
+ A = Matrix.constructWithCopy(rvals);
+ tmp = A.get(raggedr,raggedc);
+ } catch ( IllegalArgumentException e ) {
+ try_success("Catch ragged input to constructWithCopy... ",e.getMessage());
+ } catch ( java.lang.ArrayIndexOutOfBoundsException e ) {
+ errorCount = try_failure(errorCount,"Catch ragged input to constructWithCopy... ","exception not thrown in construction...ArrayIndexOutOfBoundsException thrown later");
+ }
+
+ A = new Matrix(columnwise,validld);
+ B = new Matrix(avals);
+ tmp = B.get(0,0);
+ avals[0][0] = 0.0;
+ C = B.minus(A);
+ avals[0][0] = tmp;
+ B = Matrix.constructWithCopy(avals);
+ tmp = B.get(0,0);
+ avals[0][0] = 0.0;
+ if ( ( tmp - B.get(0,0) ) != 0.0 ) {
+ /** check that constructWithCopy behaves properly **/
+ errorCount = try_failure(errorCount,"constructWithCopy... ","copy not effected... data visible outside");
+ } else {
+ try_success("constructWithCopy... ","");
+ }
+ avals[0][0] = columnwise[0];
+ I = new Matrix(ivals);
+ try {
+ check(I,Matrix.identity(3,4));
+ try_success("identity... ","");
+ } catch ( java.lang.RuntimeException e ) {
+ errorCount = try_failure(errorCount,"identity... ","identity Matrix not successfully created");
+ }
+
+/**
+ Access Methods:
+ getColumnDimension()
+ getRowDimension()
+ getArray()
+ getArrayCopy()
+ getColumnPackedCopy()
+ getRowPackedCopy()
+ get(int,int)
+ getMatrix(int,int,int,int)
+ getMatrix(int,int,int[])
+ getMatrix(int[],int,int)
+ getMatrix(int[],int[])
+ set(int,int,double)
+ setMatrix(int,int,int,int,Matrix)
+ setMatrix(int,int,int[],Matrix)
+ setMatrix(int[],int,int,Matrix)
+ setMatrix(int[],int[],Matrix)
+**/
+
+ print("\nTesting access methods...\n");
+
+/**
+ Various get methods:
+**/
+
+ B = new Matrix(avals);
+ if (B.getRowDimension() != rows) {
+ errorCount = try_failure(errorCount,"getRowDimension... ","");
+ } else {
+ try_success("getRowDimension... ","");
+ }
+ if (B.getColumnDimension() != cols) {
+ errorCount = try_failure(errorCount,"getColumnDimension... ","");
+ } else {
+ try_success("getColumnDimension... ","");
+ }
+ B = new Matrix(avals);
+ double[][] barray = B.getArray();
+ if ( barray != avals ) {
+ errorCount = try_failure(errorCount,"getArray... ","");
+ } else {
+ try_success("getArray... ","");
+ }
+ barray = B.getArrayCopy();
+ if ( barray == avals ) {
+ errorCount = try_failure(errorCount,"getArrayCopy... ","data not (deep) copied");
+ }
+ try {
+ check(barray,avals);
+ try_success("getArrayCopy... ","");
+ } catch ( java.lang.RuntimeException e ) {
+ errorCount = try_failure(errorCount,"getArrayCopy... ","data not successfully (deep) copied");
+ }
+ double[] bpacked = B.getColumnPackedCopy();
+ try {
+ check(bpacked,columnwise);
+ try_success("getColumnPackedCopy... ","");
+ } catch ( java.lang.RuntimeException e ) {
+ errorCount = try_failure(errorCount,"getColumnPackedCopy... ","data not successfully (deep) copied by columns");
+ }
+ bpacked = B.getRowPackedCopy();
+ try {
+ check(bpacked,rowwise);
+ try_success("getRowPackedCopy... ","");
+ } catch ( java.lang.RuntimeException e ) {
+ errorCount = try_failure(errorCount,"getRowPackedCopy... ","data not successfully (deep) copied by rows");
+ }
+ try {
+ tmp = B.get(B.getRowDimension(),B.getColumnDimension()-1);
+ errorCount = try_failure(errorCount,"get(int,int)... ","OutOfBoundsException expected but not thrown");
+ } catch ( java.lang.ArrayIndexOutOfBoundsException e ) {
+ try {
+ tmp = B.get(B.getRowDimension()-1,B.getColumnDimension());
+ errorCount = try_failure(errorCount,"get(int,int)... ","OutOfBoundsException expected but not thrown");
+ } catch ( java.lang.ArrayIndexOutOfBoundsException e1 ) {
+ try_success("get(int,int)... OutofBoundsException... ","");
+ }
+ } catch ( java.lang.IllegalArgumentException e1 ) {
+ errorCount = try_failure(errorCount,"get(int,int)... ","OutOfBoundsException expected but not thrown");
+ }
+ try {
+ if (B.get(B.getRowDimension()-1,B.getColumnDimension()-1) !=
+ avals[B.getRowDimension()-1][B.getColumnDimension()-1] ) {
+ errorCount = try_failure(errorCount,"get(int,int)... ","Matrix entry (i,j) not successfully retreived");
+ } else {
+ try_success("get(int,int)... ","");
+ }
+ } catch ( java.lang.ArrayIndexOutOfBoundsException e ) {
+ errorCount = try_failure(errorCount,"get(int,int)... ","Unexpected ArrayIndexOutOfBoundsException");
+ }
+ SUB = new Matrix(subavals);
+ try {
+ M = B.getMatrix(ib,ie+B.getRowDimension()+1,jb,je);
+ errorCount = try_failure(errorCount,"getMatrix(int,int,int,int)... ","ArrayIndexOutOfBoundsException expected but not thrown");
+ } catch ( java.lang.ArrayIndexOutOfBoundsException e ) {
+ try {
+ M = B.getMatrix(ib,ie,jb,je+B.getColumnDimension()+1);
+ errorCount = try_failure(errorCount,"getMatrix(int,int,int,int)... ","ArrayIndexOutOfBoundsException expected but not thrown");
+ } catch ( java.lang.ArrayIndexOutOfBoundsException e1 ) {
+ try_success("getMatrix(int,int,int,int)... ArrayIndexOutOfBoundsException... ","");
+ }
+ } catch ( java.lang.IllegalArgumentException e1 ) {
+ errorCount = try_failure(errorCount,"getMatrix(int,int,int,int)... ","ArrayIndexOutOfBoundsException expected but not thrown");
+ }
+ try {
+ M = B.getMatrix(ib,ie,jb,je);
+ try {
+ check(SUB,M);
+ try_success("getMatrix(int,int,int,int)... ","");
+ } catch ( java.lang.RuntimeException e ) {
+ errorCount = try_failure(errorCount,"getMatrix(int,int,int,int)... ","submatrix not successfully retreived");
+ }
+ } catch ( java.lang.ArrayIndexOutOfBoundsException e ) {
+ errorCount = try_failure(errorCount,"getMatrix(int,int,int,int)... ","Unexpected ArrayIndexOutOfBoundsException");
+ }
+
+ try {
+ M = B.getMatrix(ib,ie,badcolumnindexset);
+ errorCount = try_failure(errorCount,"getMatrix(int,int,int[])... ","ArrayIndexOutOfBoundsException expected but not thrown");
+ } catch ( java.lang.ArrayIndexOutOfBoundsException e ) {
+ try {
+ M = B.getMatrix(ib,ie+B.getRowDimension()+1,columnindexset);
+ errorCount = try_failure(errorCount,"getMatrix(int,int,int[])... ","ArrayIndexOutOfBoundsException expected but not thrown");
+ } catch ( java.lang.ArrayIndexOutOfBoundsException e1 ) {
+ try_success("getMatrix(int,int,int[])... ArrayIndexOutOfBoundsException... ","");
+ }
+ } catch ( java.lang.IllegalArgumentException e1 ) {
+ errorCount = try_failure(errorCount,"getMatrix(int,int,int[])... ","ArrayIndexOutOfBoundsException expected but not thrown");
+ }
+ try {
+ M = B.getMatrix(ib,ie,columnindexset);
+ try {
+ check(SUB,M);
+ try_success("getMatrix(int,int,int[])... ","");
+ } catch ( java.lang.RuntimeException e ) {
+ errorCount = try_failure(errorCount,"getMatrix(int,int,int[])... ","submatrix not successfully retreived");
+ }
+ } catch ( java.lang.ArrayIndexOutOfBoundsException e ) {
+ errorCount = try_failure(errorCount,"getMatrix(int,int,int[])... ","Unexpected ArrayIndexOutOfBoundsException");
+ }
+ try {
+ M = B.getMatrix(badrowindexset,jb,je);
+ errorCount = try_failure(errorCount,"getMatrix(int[],int,int)... ","ArrayIndexOutOfBoundsException expected but not thrown");
+ } catch ( java.lang.ArrayIndexOutOfBoundsException e ) {
+ try {
+ M = B.getMatrix(rowindexset,jb,je+B.getColumnDimension()+1);
+ errorCount = try_failure(errorCount,"getMatrix(int[],int,int)... ","ArrayIndexOutOfBoundsException expected but not thrown");
+ } catch ( java.lang.ArrayIndexOutOfBoundsException e1 ) {
+ try_success("getMatrix(int[],int,int)... ArrayIndexOutOfBoundsException... ","");
+ }
+ } catch ( java.lang.IllegalArgumentException e1 ) {
+ errorCount = try_failure(errorCount,"getMatrix(int[],int,int)... ","ArrayIndexOutOfBoundsException expected but not thrown");
+ }
+ try {
+ M = B.getMatrix(rowindexset,jb,je);
+ try {
+ check(SUB,M);
+ try_success("getMatrix(int[],int,int)... ","");
+ } catch ( java.lang.RuntimeException e ) {
+ errorCount = try_failure(errorCount,"getMatrix(int[],int,int)... ","submatrix not successfully retreived");
+ }
+ } catch ( java.lang.ArrayIndexOutOfBoundsException e ) {
+ errorCount = try_failure(errorCount,"getMatrix(int[],int,int)... ","Unexpected ArrayIndexOutOfBoundsException");
+ }
+ try {
+ M = B.getMatrix(badrowindexset,columnindexset);
+ errorCount = try_failure(errorCount,"getMatrix(int[],int[])... ","ArrayIndexOutOfBoundsException expected but not thrown");
+ } catch ( java.lang.ArrayIndexOutOfBoundsException e ) {
+ try {
+ M = B.getMatrix(rowindexset,badcolumnindexset);
+ errorCount = try_failure(errorCount,"getMatrix(int[],int[])... ","ArrayIndexOutOfBoundsException expected but not thrown");
+ } catch ( java.lang.ArrayIndexOutOfBoundsException e1 ) {
+ try_success("getMatrix(int[],int[])... ArrayIndexOutOfBoundsException... ","");
+ }
+ } catch ( java.lang.IllegalArgumentException e1 ) {
+ errorCount = try_failure(errorCount,"getMatrix(int[],int[])... ","ArrayIndexOutOfBoundsException expected but not thrown");
+ }
+ try {
+ M = B.getMatrix(rowindexset,columnindexset);
+ try {
+ check(SUB,M);
+ try_success("getMatrix(int[],int[])... ","");
+ } catch ( java.lang.RuntimeException e ) {
+ errorCount = try_failure(errorCount,"getMatrix(int[],int[])... ","submatrix not successfully retreived");
+ }
+ } catch ( java.lang.ArrayIndexOutOfBoundsException e ) {
+ errorCount = try_failure(errorCount,"getMatrix(int[],int[])... ","Unexpected ArrayIndexOutOfBoundsException");
+ }
+
+/**
+ Various set methods:
+**/
+
+ try {
+ B.set(B.getRowDimension(),B.getColumnDimension()-1,0.);
+ errorCount = try_failure(errorCount,"set(int,int,double)... ","OutOfBoundsException expected but not thrown");
+ } catch ( java.lang.ArrayIndexOutOfBoundsException e ) {
+ try {
+ B.set(B.getRowDimension()-1,B.getColumnDimension(),0.);
+ errorCount = try_failure(errorCount,"set(int,int,double)... ","OutOfBoundsException expected but not thrown");
+ } catch ( java.lang.ArrayIndexOutOfBoundsException e1 ) {
+ try_success("set(int,int,double)... OutofBoundsException... ","");
+ }
+ } catch ( java.lang.IllegalArgumentException e1 ) {
+ errorCount = try_failure(errorCount,"set(int,int,double)... ","OutOfBoundsException expected but not thrown");
+ }
+ try {
+ B.set(ib,jb,0.);
+ tmp = B.get(ib,jb);
+ try {
+ check(tmp,0.);
+ try_success("set(int,int,double)... ","");
+ } catch ( java.lang.RuntimeException e ) {
+ errorCount = try_failure(errorCount,"set(int,int,double)... ","Matrix element not successfully set");
+ }
+ } catch ( java.lang.ArrayIndexOutOfBoundsException e1) {
+ errorCount = try_failure(errorCount,"set(int,int,double)... ","Unexpected ArrayIndexOutOfBoundsException");
+ }
+ M = new Matrix(2,3,0.);
+ try {
+ B.setMatrix(ib,ie+B.getRowDimension()+1,jb,je,M);
+ errorCount = try_failure(errorCount,"setMatrix(int,int,int,int,Matrix)... ","ArrayIndexOutOfBoundsException expected but not thrown");
+ } catch ( java.lang.ArrayIndexOutOfBoundsException e ) {
+ try {
+ B.setMatrix(ib,ie,jb,je+B.getColumnDimension()+1,M);
+ errorCount = try_failure(errorCount,"setMatrix(int,int,int,int,Matrix)... ","ArrayIndexOutOfBoundsException expected but not thrown");
+ } catch ( java.lang.ArrayIndexOutOfBoundsException e1 ) {
+ try_success("setMatrix(int,int,int,int,Matrix)... ArrayIndexOutOfBoundsException... ","");
+ }
+ } catch ( java.lang.IllegalArgumentException e1 ) {
+ errorCount = try_failure(errorCount,"setMatrix(int,int,int,int,Matrix)... ","ArrayIndexOutOfBoundsException expected but not thrown");
+ }
+ try {
+ B.setMatrix(ib,ie,jb,je,M);
+ try {
+ check(M.minus(B.getMatrix(ib,ie,jb,je)),M);
+ try_success("setMatrix(int,int,int,int,Matrix)... ","");
+ } catch ( java.lang.RuntimeException e ) {
+ errorCount = try_failure(errorCount,"setMatrix(int,int,int,int,Matrix)... ","submatrix not successfully set");
+ }
+ B.setMatrix(ib,ie,jb,je,SUB);
+ } catch ( java.lang.ArrayIndexOutOfBoundsException e1 ) {
+ errorCount = try_failure(errorCount,"setMatrix(int,int,int,int,Matrix)... ","Unexpected ArrayIndexOutOfBoundsException");
+ }
+ try {
+ B.setMatrix(ib,ie+B.getRowDimension()+1,columnindexset,M);
+ errorCount = try_failure(errorCount,"setMatrix(int,int,int[],Matrix)... ","ArrayIndexOutOfBoundsException expected but not thrown");
+ } catch ( java.lang.ArrayIndexOutOfBoundsException e ) {
+ try {
+ B.setMatrix(ib,ie,badcolumnindexset,M);
+ errorCount = try_failure(errorCount,"setMatrix(int,int,int[],Matrix)... ","ArrayIndexOutOfBoundsException expected but not thrown");
+ } catch ( java.lang.ArrayIndexOutOfBoundsException e1 ) {
+ try_success("setMatrix(int,int,int[],Matrix)... ArrayIndexOutOfBoundsException... ","");
+ }
+ } catch ( java.lang.IllegalArgumentException e1 ) {
+ errorCount = try_failure(errorCount,"setMatrix(int,int,int[],Matrix)... ","ArrayIndexOutOfBoundsException expected but not thrown");
+ }
+ try {
+ B.setMatrix(ib,ie,columnindexset,M);
+ try {
+ check(M.minus(B.getMatrix(ib,ie,columnindexset)),M);
+ try_success("setMatrix(int,int,int[],Matrix)... ","");
+ } catch ( java.lang.RuntimeException e ) {
+ errorCount = try_failure(errorCount,"setMatrix(int,int,int[],Matrix)... ","submatrix not successfully set");
+ }
+ B.setMatrix(ib,ie,jb,je,SUB);
+ } catch ( java.lang.ArrayIndexOutOfBoundsException e1 ) {
+ errorCount = try_failure(errorCount,"setMatrix(int,int,int[],Matrix)... ","Unexpected ArrayIndexOutOfBoundsException");
+ }
+ try {
+ B.setMatrix(rowindexset,jb,je+B.getColumnDimension()+1,M);
+ errorCount = try_failure(errorCount,"setMatrix(int[],int,int,Matrix)... ","ArrayIndexOutOfBoundsException expected but not thrown");
+ } catch ( java.lang.ArrayIndexOutOfBoundsException e ) {
+ try {
+ B.setMatrix(badrowindexset,jb,je,M);
+ errorCount = try_failure(errorCount,"setMatrix(int[],int,int,Matrix)... ","ArrayIndexOutOfBoundsException expected but not thrown");
+ } catch ( java.lang.ArrayIndexOutOfBoundsException e1 ) {
+ try_success("setMatrix(int[],int,int,Matrix)... ArrayIndexOutOfBoundsException... ","");
+ }
+ } catch ( java.lang.IllegalArgumentException e1 ) {
+ errorCount = try_failure(errorCount,"setMatrix(int[],int,int,Matrix)... ","ArrayIndexOutOfBoundsException expected but not thrown");
+ }
+ try {
+ B.setMatrix(rowindexset,jb,je,M);
+ try {
+ check(M.minus(B.getMatrix(rowindexset,jb,je)),M);
+ try_success("setMatrix(int[],int,int,Matrix)... ","");
+ } catch ( java.lang.RuntimeException e ) {
+ errorCount = try_failure(errorCount,"setMatrix(int[],int,int,Matrix)... ","submatrix not successfully set");
+ }
+ B.setMatrix(ib,ie,jb,je,SUB);
+ } catch ( java.lang.ArrayIndexOutOfBoundsException e1 ) {
+ errorCount = try_failure(errorCount,"setMatrix(int[],int,int,Matrix)... ","Unexpected ArrayIndexOutOfBoundsException");
+ }
+ try {
+ B.setMatrix(rowindexset,badcolumnindexset,M);
+ errorCount = try_failure(errorCount,"setMatrix(int[],int[],Matrix)... ","ArrayIndexOutOfBoundsException expected but not thrown");
+ } catch ( java.lang.ArrayIndexOutOfBoundsException e ) {
+ try {
+ B.setMatrix(badrowindexset,columnindexset,M);
+ errorCount = try_failure(errorCount,"setMatrix(int[],int[],Matrix)... ","ArrayIndexOutOfBoundsException expected but not thrown");
+ } catch ( java.lang.ArrayIndexOutOfBoundsException e1 ) {
+ try_success("setMatrix(int[],int[],Matrix)... ArrayIndexOutOfBoundsException... ","");
+ }
+ } catch ( java.lang.IllegalArgumentException e1 ) {
+ errorCount = try_failure(errorCount,"setMatrix(int[],int[],Matrix)... ","ArrayIndexOutOfBoundsException expected but not thrown");
+ }
+ try {
+ B.setMatrix(rowindexset,columnindexset,M);
+ try {
+ check(M.minus(B.getMatrix(rowindexset,columnindexset)),M);
+ try_success("setMatrix(int[],int[],Matrix)... ","");
+ } catch ( java.lang.RuntimeException e ) {
+ errorCount = try_failure(errorCount,"setMatrix(int[],int[],Matrix)... ","submatrix not successfully set");
+ }
+ } catch ( java.lang.ArrayIndexOutOfBoundsException e1 ) {
+ errorCount = try_failure(errorCount,"setMatrix(int[],int[],Matrix)... ","Unexpected ArrayIndexOutOfBoundsException");
+ }
+
+/**
+ Array-like methods:
+ minus
+ minusEquals
+ plus
+ plusEquals
+ arrayLeftDivide
+ arrayLeftDivideEquals
+ arrayRightDivide
+ arrayRightDivideEquals
+ arrayTimes
+ arrayTimesEquals
+ uminus
+**/
+
+ print("\nTesting array-like methods...\n");
+ S = new Matrix(columnwise,nonconformld);
+ R = Matrix.random(A.getRowDimension(),A.getColumnDimension());
+ A = R;
+ try {
+ S = A.minus(S);
+ errorCount = try_failure(errorCount,"minus conformance check... ","nonconformance not raised");
+ } catch ( IllegalArgumentException e ) {
+ try_success("minus conformance check... ","");
+ }
+ if (A.minus(R).norm1() != 0.) {
+ errorCount = try_failure(errorCount,"minus... ","(difference of identical Matrices is nonzero,\nSubsequent use of minus should be suspect)");
+ } else {
+ try_success("minus... ","");
+ }
+ A = R.copy();
+ A.minusEquals(R);
+ Z = new Matrix(A.getRowDimension(),A.getColumnDimension());
+ try {
+ A.minusEquals(S);
+ errorCount = try_failure(errorCount,"minusEquals conformance check... ","nonconformance not raised");
+ } catch ( IllegalArgumentException e ) {
+ try_success("minusEquals conformance check... ","");
+ }
+ if (A.minus(Z).norm1() != 0.) {
+ errorCount = try_failure(errorCount,"minusEquals... ","(difference of identical Matrices is nonzero,\nSubsequent use of minus should be suspect)");
+ } else {
+ try_success("minusEquals... ","");
+ }
+
+ A = R.copy();
+ B = Matrix.random(A.getRowDimension(),A.getColumnDimension());
+ C = A.minus(B);
+ try {
+ S = A.plus(S);
+ errorCount = try_failure(errorCount,"plus conformance check... ","nonconformance not raised");
+ } catch ( IllegalArgumentException e ) {
+ try_success("plus conformance check... ","");
+ }
+ try {
+ check(C.plus(B),A);
+ try_success("plus... ","");
+ } catch ( java.lang.RuntimeException e ) {
+ errorCount = try_failure(errorCount,"plus... ","(C = A - B, but C + B != A)");
+ }
+ C = A.minus(B);
+ C.plusEquals(B);
+ try {
+ A.plusEquals(S);
+ errorCount = try_failure(errorCount,"plusEquals conformance check... ","nonconformance not raised");
+ } catch ( IllegalArgumentException e ) {
+ try_success("plusEquals conformance check... ","");
+ }
+ try {
+ check(C,A);
+ try_success("plusEquals... ","");
+ } catch ( java.lang.RuntimeException e ) {
+ errorCount = try_failure(errorCount,"plusEquals... ","(C = A - B, but C = C + B != A)");
+ }
+ A = R.uminus();
+ try {
+ check(A.plus(R),Z);
+ try_success("uminus... ","");
+ } catch ( java.lang.RuntimeException e ) {
+ errorCount = try_failure(errorCount,"uminus... ","(-A + A != zeros)");
+ }
+ A = R.copy();
+ O = new Matrix(A.getRowDimension(),A.getColumnDimension(),1.0);
+ C = A.arrayLeftDivide(R);
+ try {
+ S = A.arrayLeftDivide(S);
+ errorCount = try_failure(errorCount,"arrayLeftDivide conformance check... ","nonconformance not raised");
+ } catch ( IllegalArgumentException e ) {
+ try_success("arrayLeftDivide conformance check... ","");
+ }
+ try {
+ check(C,O);
+ try_success("arrayLeftDivide... ","");
+ } catch ( java.lang.RuntimeException e ) {
+ errorCount = try_failure(errorCount,"arrayLeftDivide... ","(M.\\M != ones)");
+ }
+ try {
+ A.arrayLeftDivideEquals(S);
+ errorCount = try_failure(errorCount,"arrayLeftDivideEquals conformance check... ","nonconformance not raised");
+ } catch ( IllegalArgumentException e ) {
+ try_success("arrayLeftDivideEquals conformance check... ","");
+ }
+ A.arrayLeftDivideEquals(R);
+ try {
+ check(A,O);
+ try_success("arrayLeftDivideEquals... ","");
+ } catch ( java.lang.RuntimeException e ) {
+ errorCount = try_failure(errorCount,"arrayLeftDivideEquals... ","(M.\\M != ones)");
+ }
+ A = R.copy();
+ try {
+ A.arrayRightDivide(S);
+ errorCount = try_failure(errorCount,"arrayRightDivide conformance check... ","nonconformance not raised");
+ } catch ( IllegalArgumentException e ) {
+ try_success("arrayRightDivide conformance check... ","");
+ }
+ C = A.arrayRightDivide(R);
+ try {
+ check(C,O);
+ try_success("arrayRightDivide... ","");
+ } catch ( java.lang.RuntimeException e ) {
+ errorCount = try_failure(errorCount,"arrayRightDivide... ","(M./M != ones)");
+ }
+ try {
+ A.arrayRightDivideEquals(S);
+ errorCount = try_failure(errorCount,"arrayRightDivideEquals conformance check... ","nonconformance not raised");
+ } catch ( IllegalArgumentException e ) {
+ try_success("arrayRightDivideEquals conformance check... ","");
+ }
+ A.arrayRightDivideEquals(R);
+ try {
+ check(A,O);
+ try_success("arrayRightDivideEquals... ","");
+ } catch ( java.lang.RuntimeException e ) {
+ errorCount = try_failure(errorCount,"arrayRightDivideEquals... ","(M./M != ones)");
+ }
+ A = R.copy();
+ B = Matrix.random(A.getRowDimension(),A.getColumnDimension());
+ try {
+ S = A.arrayTimes(S);
+ errorCount = try_failure(errorCount,"arrayTimes conformance check... ","nonconformance not raised");
+ } catch ( IllegalArgumentException e ) {
+ try_success("arrayTimes conformance check... ","");
+ }
+ C = A.arrayTimes(B);
+ try {
+ check(C.arrayRightDivideEquals(B),A);
+ try_success("arrayTimes... ","");
+ } catch ( java.lang.RuntimeException e ) {
+ errorCount = try_failure(errorCount,"arrayTimes... ","(A = R, C = A.*B, but C./B != A)");
+ }
+ try {
+ A.arrayTimesEquals(S);
+ errorCount = try_failure(errorCount,"arrayTimesEquals conformance check... ","nonconformance not raised");
+ } catch ( IllegalArgumentException e ) {
+ try_success("arrayTimesEquals conformance check... ","");
+ }
+ A.arrayTimesEquals(B);
+ try {
+ check(A.arrayRightDivideEquals(B),R);
+ try_success("arrayTimesEquals... ","");
+ } catch ( java.lang.RuntimeException e ) {
+ errorCount = try_failure(errorCount,"arrayTimesEquals... ","(A = R, A = A.*B, but A./B != R)");
+ }
+
+/**
+ I/O methods:
+ read
+ print
+ serializable:
+ writeObject
+ readObject
+**/
+ print("\nTesting I/O methods...\n");
+ try {
+ DecimalFormat fmt = new DecimalFormat("0.0000E00");
+ fmt.setDecimalFormatSymbols(new DecimalFormatSymbols(Locale.US));
+
+ PrintWriter FILE = new PrintWriter(new FileOutputStream("JamaTestMatrix.out"));
+ A.print(FILE,fmt,10);
+ FILE.close();
+ R = Matrix.read(new BufferedReader(new FileReader("JamaTestMatrix.out")));
+ if (A.minus(R).norm1() < .001 ) {
+ try_success("print()/read()...","");
+ } else {
+ errorCount = try_failure(errorCount,"print()/read()...","Matrix read from file does not match Matrix printed to file");
+ }
+ } catch ( java.io.IOException ioe ) {
+ warningCount = try_warning(warningCount,"print()/read()...","unexpected I/O error, unable to run print/read test; check write permission in current directory and retry");
+ } catch(Exception e) {
+ try {
+ e.printStackTrace(System.out);
+ warningCount = try_warning(warningCount,"print()/read()...","Formatting error... will try JDK1.1 reformulation...");
+ DecimalFormat fmt = new DecimalFormat("0.0000");
+ PrintWriter FILE = new PrintWriter(new FileOutputStream("JamaTestMatrix.out"));
+ A.print(FILE,fmt,10);
+ FILE.close();
+ R = Matrix.read(new BufferedReader(new FileReader("JamaTestMatrix.out")));
+ if (A.minus(R).norm1() < .001 ) {
+ try_success("print()/read()...","");
+ } else {
+ errorCount = try_failure(errorCount,"print()/read() (2nd attempt) ...","Matrix read from file does not match Matrix printed to file");
+ }
+ } catch ( java.io.IOException ioe ) {
+ warningCount = try_warning(warningCount,"print()/read()...","unexpected I/O error, unable to run print/read test; check write permission in current directory and retry");
+ }
+ }
+
+ R = Matrix.random(A.getRowDimension(),A.getColumnDimension());
+ String tmpname = "TMPMATRIX.serial";
+ try {
+ ObjectOutputStream out = new ObjectOutputStream(new FileOutputStream(tmpname));
+ out.writeObject(R);
+ ObjectInputStream sin = new ObjectInputStream(new FileInputStream(tmpname));
+ A = (Matrix) sin.readObject();
+
+ try {
+ check(A,R);
+ try_success("writeObject(Matrix)/readObject(Matrix)...","");
+ } catch ( java.lang.RuntimeException e ) {
+ errorCount = try_failure(errorCount,"writeObject(Matrix)/readObject(Matrix)...","Matrix not serialized correctly");
+ }
+ } catch ( java.io.IOException ioe ) {
+ warningCount = try_warning(warningCount,"writeObject()/readObject()...","unexpected I/O error, unable to run serialization test; check write permission in current directory and retry");
+ } catch(Exception e) {
+ errorCount = try_failure(errorCount,"writeObject(Matrix)/readObject(Matrix)...","unexpected error in serialization test");
+ }
+
+/**
+ LA methods:
+ transpose
+ times
+ cond
+ rank
+ det
+ trace
+ norm1
+ norm2
+ normF
+ normInf
+ solve
+ solveTranspose
+ inverse
+ chol
+ eig
+ lu
+ qr
+ svd
+**/
+
+ print("\nTesting linear algebra methods...\n");
+ A = new Matrix(columnwise,3);
+ T = new Matrix(tvals);
+ T = A.transpose();
+ try {
+ check(A.transpose(),T);
+ try_success("transpose...","");
+ } catch ( java.lang.RuntimeException e ) {
+ errorCount = try_failure(errorCount,"transpose()...","transpose unsuccessful");
+ }
+ A.transpose();
+ try {
+ check(A.norm1(),columnsummax);
+ try_success("norm1...","");
+ } catch ( java.lang.RuntimeException e ) {
+ errorCount = try_failure(errorCount,"norm1()...","incorrect norm calculation");
+ }
+ try {
+ check(A.normInf(),rowsummax);
+ try_success("normInf()...","");
+ } catch ( java.lang.RuntimeException e ) {
+ errorCount = try_failure(errorCount,"normInf()...","incorrect norm calculation");
+ }
+ try {
+ check(A.normF(),Math.sqrt(sumofsquares));
+ try_success("normF...","");
+ } catch ( java.lang.RuntimeException e ) {
+ errorCount = try_failure(errorCount,"normF()...","incorrect norm calculation");
+ }
+ try {
+ check(A.trace(),sumofdiagonals);
+ try_success("trace()...","");
+ } catch ( java.lang.RuntimeException e ) {
+ errorCount = try_failure(errorCount,"trace()...","incorrect trace calculation");
+ }
+ try {
+ check(A.getMatrix(0,A.getRowDimension()-1,0,A.getRowDimension()-1).det(),0.);
+ try_success("det()...","");
+ } catch ( java.lang.RuntimeException e ) {
+ errorCount = try_failure(errorCount,"det()...","incorrect determinant calculation");
+ }
+ SQ = new Matrix(square);
+ try {
+ check(A.times(A.transpose()),SQ);
+ try_success("times(Matrix)...","");
+ } catch ( java.lang.RuntimeException e ) {
+ errorCount = try_failure(errorCount,"times(Matrix)...","incorrect Matrix-Matrix product calculation");
+ }
+ try {
+ check(A.times(0.),Z);
+ try_success("times(double)...","");
+ } catch ( java.lang.RuntimeException e ) {
+ errorCount = try_failure(errorCount,"times(double)...","incorrect Matrix-scalar product calculation");
+ }
+
+ A = new Matrix(columnwise,4);
+ QRDecomposition QR = A.qr();
+ R = QR.getR();
+ try {
+ check(A,QR.getQ().times(R));
+ try_success("QRDecomposition...","");
+ } catch ( java.lang.RuntimeException e ) {
+ errorCount = try_failure(errorCount,"QRDecomposition...","incorrect QR decomposition calculation");
+ }
+ SingularValueDecomposition SVD = A.svd();
+ try {
+ check(A,SVD.getU().times(SVD.getS().times(SVD.getV().transpose())));
+ try_success("SingularValueDecomposition...","");
+ } catch ( java.lang.RuntimeException e ) {
+ errorCount = try_failure(errorCount,"SingularValueDecomposition...","incorrect singular value decomposition calculation");
+ }
+ DEF = new Matrix(rankdef);
+ try {
+ check(DEF.rank(),Math.min(DEF.getRowDimension(),DEF.getColumnDimension())-1);
+ try_success("rank()...","");
+ } catch ( java.lang.RuntimeException e ) {
+ errorCount = try_failure(errorCount,"rank()...","incorrect rank calculation");
+ }
+ B = new Matrix(condmat);
+ SVD = B.svd();
+ double [] singularvalues = SVD.getSingularValues();
+ try {
+ check(B.cond(),singularvalues[0]/singularvalues[Math.min(B.getRowDimension(),B.getColumnDimension())-1]);
+ try_success("cond()...","");
+ } catch ( java.lang.RuntimeException e ) {
+ errorCount = try_failure(errorCount,"cond()...","incorrect condition number calculation");
+ }
+ int n = A.getColumnDimension();
+ A = A.getMatrix(0,n-1,0,n-1);
+ A.set(0,0,0.);
+ LUDecomposition LU = A.lu();
+ try {
+ check(A.getMatrix(LU.getPivot(),0,n-1),LU.getL().times(LU.getU()));
+ try_success("LUDecomposition...","");
+ } catch ( java.lang.RuntimeException e ) {
+ errorCount = try_failure(errorCount,"LUDecomposition...","incorrect LU decomposition calculation");
+ }
+ X = A.inverse();
+ try {
+ check(A.times(X),Matrix.identity(3,3));
+ try_success("inverse()...","");
+ } catch ( java.lang.RuntimeException e ) {
+ errorCount = try_failure(errorCount,"inverse()...","incorrect inverse calculation");
+ }
+ O = new Matrix(SUB.getRowDimension(),1,1.0);
+ SOL = new Matrix(sqSolution);
+ SQ = SUB.getMatrix(0,SUB.getRowDimension()-1,0,SUB.getRowDimension()-1);
+ try {
+ check(SQ.solve(SOL),O);
+ try_success("solve()...","");
+ } catch ( java.lang.IllegalArgumentException e1 ) {
+ errorCount = try_failure(errorCount,"solve()...",e1.getMessage());
+ } catch ( java.lang.RuntimeException e ) {
+ errorCount = try_failure(errorCount,"solve()...",e.getMessage());
+ }
+ A = new Matrix(pvals);
+ CholeskyDecomposition Chol = A.chol();
+ Matrix L = Chol.getL();
+ try {
+ check(A,L.times(L.transpose()));
+ try_success("CholeskyDecomposition...","");
+ } catch ( java.lang.RuntimeException e ) {
+ errorCount = try_failure(errorCount,"CholeskyDecomposition...","incorrect Cholesky decomposition calculation");
+ }
+ X = Chol.solve(Matrix.identity(3,3));
+ try {
+ check(A.times(X),Matrix.identity(3,3));
+ try_success("CholeskyDecomposition solve()...","");
+ } catch ( java.lang.RuntimeException e ) {
+ errorCount = try_failure(errorCount,"CholeskyDecomposition solve()...","incorrect Choleskydecomposition solve calculation");
+ }
+ EigenvalueDecomposition Eig = A.eig();
+ Matrix D = Eig.getD();
+ Matrix V = Eig.getV();
+ try {
+ check(A.times(V),V.times(D));
+ try_success("EigenvalueDecomposition (symmetric)...","");
+ } catch ( java.lang.RuntimeException e ) {
+ errorCount = try_failure(errorCount,"EigenvalueDecomposition (symmetric)...","incorrect symmetric Eigenvalue decomposition calculation");
+ }
+ A = new Matrix(evals);
+ Eig = A.eig();
+ D = Eig.getD();
+ V = Eig.getV();
+ try {
+ check(A.times(V),V.times(D));
+ try_success("EigenvalueDecomposition (nonsymmetric)...","");
+ } catch ( java.lang.RuntimeException e ) {
+ errorCount = try_failure(errorCount,"EigenvalueDecomposition (nonsymmetric)...","incorrect nonsymmetric Eigenvalue decomposition calculation");
+ }
+
+ try {
+ print("\nTesting Eigenvalue; If this hangs, we've failed\n");
+ Matrix bA = new Matrix(badeigs);
+ EigenvalueDecomposition bEig = bA.eig();
+ try_success("EigenvalueDecomposition (hang)...","");
+ } catch ( java.lang.RuntimeException e ) {
+ errorCount = try_failure(errorCount,"EigenvalueDecomposition (hang)...",
+ "incorrect termination");
+ }
+
+
+ print("\nTestMatrix completed.\n");
+ print("Total errors reported: " + Integer.toString(errorCount) + "\n");
+ print("Total warnings reported: " + Integer.toString(warningCount) + "\n");
+ }
+
+ /** private utility routines **/
+
+ /** Check magnitude of difference of scalars. **/
+
+ private static void check(double x, double y) {
+ double eps = Math.pow(2.0,-52.0);
+ if (x == 0 & Math.abs(y) < 10*eps) return;
+ if (y == 0 & Math.abs(x) < 10*eps) return;
+ if (Math.abs(x-y) > 10*eps*Math.max(Math.abs(x),Math.abs(y))) {
+ throw new RuntimeException("The difference x-y is too large: x = " + Double.toString(x) + " y = " + Double.toString(y));
+ }
+ }
+
+ /** Check norm of difference of "vectors". **/
+
+ private static void check(double[] x, double[] y) {
+ if (x.length == y.length ) {
+ for (int i=0;i 1000*eps*Math.max(X.norm1(),Y.norm1())) {
+ throw new RuntimeException("The norm of (X-Y) is too large: " + Double.toString(X.minus(Y).norm1()));
+ }
+ }
+
+ /** Shorten spelling of print. **/
+
+ private static void print (String s) {
+ System.out.print(s);
+ }
+
+ /** Print appropriate messages for successful outcome try **/
+
+ private static void try_success (String s,String e) {
+ print("> " + s + "success\n");
+ if ( e != "" ) {
+ print("> Message: " + e + "\n");
+ }
+ }
+ /** Print appropriate messages for unsuccessful outcome try **/
+
+ private static int try_failure (int count, String s,String e) {
+ print("> " + s + "*** failure ***\n> Message: " + e + "\n");
+ return ++count;
+ }
+
+ /** Print appropriate messages for unsuccessful outcome try **/
+
+ private static int try_warning (int count, String s,String e) {
+ print("> " + s + "*** warning ***\n> Message: " + e + "\n");
+ return ++count;
+ }
+
+ /** Print a row vector. **/
+
+ private static void print(double[] x, int w, int d) {
+ // Use format Fw.d for all elements.
+ System.out.print("\n");
+ new Matrix(x,1).print(w,d);
+ print("\n");
+ }
+
+}
diff --git a/Java/src/Jama/util/Maths.java b/Java/src/Jama/util/Maths.java
new file mode 100644
index 00000000..087a6a74
--- /dev/null
+++ b/Java/src/Jama/util/Maths.java
@@ -0,0 +1,20 @@
+package Jama.util;
+
+public class Maths {
+
+ /** sqrt(a^2 + b^2) without under/overflow. **/
+
+ public static double hypot(double a, double b) {
+ double r;
+ if (Math.abs(a) > Math.abs(b)) {
+ r = b/a;
+ r = Math.abs(a)*Math.sqrt(1+r*r);
+ } else if (b != 0) {
+ r = a/b;
+ r = Math.abs(b)*Math.sqrt(1+r*r);
+ } else {
+ r = 0.0;
+ }
+ return r;
+ }
+}
diff --git a/Java/src/Viewer2D/EdgesFrame.java b/Java/src/Viewer2D/EdgesFrame.java
index b73c609c..0cde000c 100755
--- a/Java/src/Viewer2D/EdgesFrame.java
+++ b/Java/src/Viewer2D/EdgesFrame.java
@@ -21,22 +21,32 @@
import java.awt.BorderLayout;
import java.awt.Color;
import java.awt.Container;
+import java.awt.Cursor;
import java.awt.Dimension;
import java.awt.FlowLayout;
import java.awt.Font;
import java.awt.Graphics;
+import java.awt.Image;
+import java.awt.Point;
import java.awt.Toolkit;
import java.awt.event.ActionEvent;
import java.awt.event.ActionListener;
import java.awt.event.ItemEvent;
import java.awt.event.ItemListener;
+import java.awt.event.KeyListener;
+import java.awt.event.ComponentListener;
+import java.awt.event.ComponentEvent;
+import java.awt.event.KeyEvent;
+import java.awt.event.InputEvent;
import java.io.File;
import java.io.FileNotFoundException;
import java.io.IOException;
import java.util.HashMap;
+import javax.swing.AbstractAction;
import javax.swing.JButton;
import javax.swing.JColorChooser;
+import javax.swing.JComponent;
import javax.swing.JFileChooser;
import javax.swing.JFrame;
import javax.swing.JMenu;
@@ -44,10 +54,15 @@
import javax.swing.JMenuItem;
import javax.swing.JOptionPane;
import javax.swing.JPanel;
+import javax.swing.JPopupMenu;
import javax.swing.JRadioButton;
import javax.swing.JRadioButtonMenuItem;
import javax.swing.JScrollPane;
import javax.swing.JTextField;
+import javax.swing.JToggleButton;
+import javax.swing.KeyStroke;
+
+import Viewer2D.EdgesPanel.modes;
/**
* SESS - 2014.05.11:
@@ -66,7 +81,7 @@
// TODO: SESS - Also show "working, please wait..."
// TODO: SESS - When edges panel should repaint is repainting this frame
// should be only the edges panel
-public class EdgesFrame extends JFrame {
+public class EdgesFrame extends JFrame implements KeyListener,ComponentListener {
private static final long serialVersionUID = 1099866981814538683L;
static final int FRAME_WIDTH = 1024;
@@ -75,6 +90,10 @@ public class EdgesFrame extends JFrame {
private Edge[] edges;
private Vertex[] vertices;
+ private HashMap labels;
+
+ private double labelScale;
+
private EdgesPanel panel; // Has the edges drawn
private JTextField statusBar; // Shows highlighted/user info
@@ -119,6 +138,11 @@ public class EdgesFrame extends JFrame {
private Color fontColor, vertexColor, edgeColor, backgroundColor;
private int STATUSBAR_Y;
+
+
+ JToggleButton hand;
+
+ JToggleButton magnifier;
public EdgesFrame(String title, int x, int y) {
super(title);
@@ -128,7 +152,7 @@ public EdgesFrame(String title, int x, int y) {
STATUSBAR_Y = 50;
moveStepSize = .2;
- zoomStepSize = .2;
+ zoomStepSize = .8;
windowSizes = new int[2];
windowSizes[0] = x;
@@ -147,12 +171,16 @@ public EdgesFrame(String title, int x, int y) {
container.setLayout(layout);
setDefaultCloseOperation(JFrame.EXIT_ON_CLOSE);
-
- panel = new EdgesPanel(edges, vertices, windowSizes[0], windowSizes[1]);
+ Dimension screen = Toolkit.getDefaultToolkit().getScreenSize();
+ panel = new EdgesPanel(edges, vertices, new HashMap(),screen.width,screen.height);//windowSizes[0], windowSizes[1]);
panel.setBackground(backgroundColor);
setMenuBars();
+ hand = new JToggleButton("Hand");
+ magnifier = new JToggleButton("Magnifier");
setButtons();
+
+ addKeyListener(this);
JScrollPane scrollPane = new JScrollPane(panel);
container.add(scrollPane, BorderLayout.CENTER);
@@ -163,13 +191,85 @@ public EdgesFrame(String title, int x, int y) {
panel.setStatusBar(statusBar);
// Let's start with the screen size
- Dimension screen = Toolkit.getDefaultToolkit().getScreenSize();
+
int width = Math.min(screen.width, FRAME_WIDTH);
int height = Math.min(screen.height, FRAME_HEIGHT);
setSize(width, height);
setLocation((screen.width - width) / 2, (screen.height - height) / 2);
setResizable(true);
setVisible(true);
+ labels = new HashMap<>();
+ labelScale = 1;
+ addComponentListener(this);
+ /*ActionListener okListener = new ActionListener() {
+
+ @Override
+ public void actionPerformed(ActionEvent e) {
+ System.out.println("2");
+ }
+ };
+ getRootPane().registerKeyboardAction(okListener,
+ KeyStroke.getKeyStroke(KeyEvent.VK_SPACE, 0),
+ JComponent.WHEN_IN_FOCUSED_WINDOW);*/
+ //KeyStroke k = .getKeyStroke(' ' );
+ //k.getKeyStroke(keyChar, onKeyRelease)
+
+
+ panel.getInputMap(JComponent.WHEN_IN_FOCUSED_WINDOW).put(KeyStroke.getKeyStroke(' ' ), "SPACE");
+
+ panel.getActionMap().put("SPACE", new AbstractAction() {
+ public void actionPerformed(ActionEvent e) {
+ //do nothing
+
+ if (panel.getMode()!=modes.handmode)
+ {
+ setMode(modes.handmode);
+ }
+ else
+ setMode(modes.nomode);
+ }
+ });
+
+
+ panel.getInputMap(JComponent.WHEN_IN_FOCUSED_WINDOW).put(KeyStroke.getKeyStroke( KeyEvent.VK_ADD, InputEvent.CTRL_MASK ), "ZOOMIN");
+
+ panel.getInputMap(JComponent.WHEN_IN_FOCUSED_WINDOW).put(KeyStroke.getKeyStroke(KeyEvent.VK_EQUALS, KeyEvent.CTRL_DOWN_MASK | KeyEvent.SHIFT_DOWN_MASK),"ZOOMIN");
+
+
+ panel.getActionMap().put("ZOOMIN", new AbstractAction() {
+ public void actionPerformed(ActionEvent e) {
+ Point p = panel.getMousePosition();
+ int x1 = (int)p.getX();
+ int y1 = (int)p.getY();
+
+ panel.zoom2Point(x1, y1);
+ panel.setPaintImage();
+ panel.repaint();
+ }
+ });
+
+
+
+ panel.getInputMap(JComponent.WHEN_IN_FOCUSED_WINDOW).put(KeyStroke.getKeyStroke(KeyEvent.VK_SUBTRACT, InputEvent.CTRL_MASK ), "ZOOMOUT");
+
+ //panel.getInputMap(JComponent.WHEN_IN_FOCUSED_WINDOW).put(KeyStroke.getKeyStroke(KeyEvent.VK_MINUS, KeyEvent.CTRL_DOWN_MASK | KeyEvent.SHIFT_DOWN_MASK),"ZOOMOUT");
+
+ panel.getInputMap(JComponent.WHEN_IN_FOCUSED_WINDOW).put(KeyStroke.getKeyStroke(KeyEvent.VK_MINUS, KeyEvent.CTRL_DOWN_MASK),"ZOOMOUT");
+
+
+ panel.getActionMap().put("ZOOMOUT", new AbstractAction() {
+ public void actionPerformed(ActionEvent e) {
+ Point p = panel.getMousePosition();
+ int x1 = (int)p.getX();
+ int y1 = (int)p.getY();
+
+ panel.zoomOutFromPoint(x1, y1);
+ panel.setPaintImage();
+ panel.repaint();
+ }
+ });
+
+
}
private JFileChooser fileChooser;
@@ -210,9 +310,11 @@ private void setMenuBars() {
JMenu fileMenu = new JMenu("File");
fileMenu.setMnemonic('F');
- // LOAD THE EDGES FILE ( SOMETHING.lgl )
- JMenuItem openFile = new JMenuItem("Open .lgl file");
- openFile.setMnemonic('O');
+
+ // load all
+
+ JMenuItem openFile = new JMenuItem("Load all files");
+
openFile.addActionListener(new ActionListener() {
public void actionPerformed(ActionEvent e) {
JFileChooser chooser = getFileChooser("lgl");
@@ -220,11 +322,38 @@ public void actionPerformed(ActionEvent e) {
if (returnVal == JFileChooser.APPROVE_OPTION) {
edgesFile = chooser.getSelectedFile();
loadSHORTFile(edgesFile);
+ File coordsFile = new File(edgesFile.getPath().replace(".lgl", ".coords"));
+ if (coordsFile.isFile())
+ loadCoordsFile(coordsFile);
+ File colorsFile = new File(edgesFile.getPath().replace(".lgl", ".colors"));
+
+ if (colorsFile.isFile())
+ loadEdgeColorFile(colorsFile);
+
+ File labelsFile = new File(edgesFile.getPath().replace(".lgl", ".labels"));
+
+ if (labelsFile.isFile())
+ loadLabelFile(labelsFile);
}
}
});
fileMenu.add(openFile);
+ // LOAD THE EDGES FILE ( SOMETHING.lgl )
+ JMenuItem openAllFile = new JMenuItem("Open .lgl file");
+ openAllFile.setMnemonic('O');
+ openAllFile.addActionListener(new ActionListener() {
+ public void actionPerformed(ActionEvent e) {
+ JFileChooser chooser = getFileChooser("lgl");
+ int returnVal = chooser.showOpenDialog(EdgesFrame.this);
+ if (returnVal == JFileChooser.APPROVE_OPTION) {
+ edgesFile = chooser.getSelectedFile();
+ loadSHORTFile(edgesFile);
+ }
+ }
+ });
+ fileMenu.add(openAllFile);
+
// LOAD THE 2D COORDS
JMenuItem cFile = new JMenuItem("Open 2D Coords file");
cFile.setMnemonic('C');
@@ -335,6 +464,27 @@ public void actionPerformed(ActionEvent e) {
});
fileMenu.add(reload);
+ // load labels
+
+ JMenuItem labelMenu = new JMenuItem("Load labels");
+ labelMenu.setMnemonic('L');
+ labelMenu.addActionListener(new ActionListener() {
+ public void actionPerformed(ActionEvent e) {
+ if (vertexCheck()) {
+ JFileChooser chooser = getFileChooser();
+ int returnVal = chooser.showOpenDialog(EdgesFrame.this);
+ if (returnVal == JFileChooser.APPROVE_OPTION) {
+ loadLabelFile(chooser.getSelectedFile());
+ }
+ }
+ }
+ });
+ fileMenu.add(labelMenu);
+
+
+
+
+
// APPLICATION EXIT
JMenuItem exit = new JMenuItem("Exit");
exit.setMnemonic('x');
@@ -391,13 +541,42 @@ public void actionPerformed(ActionEvent e) {
"Error", JOptionPane.ERROR_MESSAGE);
} else {
zoomStepSize = newZoomStepSize;
- panel.setZoomStepSize(zoomStepSize);
+ formatter.setLabelScale(newZoomStepSize);
+
+
}
}
}
});
edit.add(zoomsize);
+ // change label scale
+
+ JMenuItem labelscalesize = new JMenuItem("Change scale for Labels");
+
+ labelscalesize.addActionListener(new ActionListener() {
+ public void actionPerformed(ActionEvent e) {
+ String newSize = JOptionPane.showInputDialog(EdgesFrame.this,
+ "Enter the new label scale",
+ Double.toString(labelScale));
+ if (newSize != null) {
+ double newLabelScale = Double.parseDouble(newSize);
+ if (newLabelScale <= 0 ) {
+ JOptionPane.showMessageDialog(null, "Illegal Value",
+ "Error", JOptionPane.ERROR_MESSAGE);
+ } else {
+ labelScale = newLabelScale;
+
+ formatter.setLabelScale(labelScale);
+ panel.setPaintImage(); // TODO: move inside panel.xxxx()?
+ panel.repaint();
+ }
+ }
+ }
+ });
+ edit.add(labelscalesize);
+
+
// BUTTON TO REMOVE TRANSIENT EDGES
blockers = new JRadioButton("Remove Transient Edges");
blockers.setMnemonic('R');
@@ -433,6 +612,8 @@ public void actionPerformed(ActionEvent e) {
if (zoomHighlighted.isSelected()) {
panel.prepZoomPoint(true);
statusMessage = "Enabled Zoom Point";
+ setMode(modes.nomode);
+
} else {
panel.prepZoomPoint(false);
statusMessage = "Disabled Zoom Point";
@@ -451,6 +632,7 @@ public void actionPerformed(ActionEvent e) {
if (zoomRegion.isSelected()) {
panel.prepZoomRegion(true);
statusMessage = "Enabled Zoom Region";
+ setMode(modes.nomode);
} else {
panel.prepZoomRegion(false);
statusMessage = "Disabled Zoom Region";
@@ -802,6 +984,23 @@ public void loadVertexColorFile(File f) {
}
}
+ public void loadLabelFile(File f) {
+ System.out.println("Loading Label File " + f.getAbsolutePath());
+ try {
+ edgesio.loadLabelFile(f);
+ // Update the panel with tlabels
+ panel.addLabels(edgesio.getLabels());
+ panel.setPaintImage();
+ panel.repaint();
+ } catch (FileNotFoundException ee) {
+ JOptionPane.showMessageDialog(null, "File Not Found", "Error",
+ JOptionPane.ERROR_MESSAGE);
+ } catch (IOException ee) {
+ JOptionPane.showMessageDialog(null, "IO Error: Check File Format",
+ "Error", JOptionPane.ERROR_MESSAGE);
+ }
+ }
+
public void loadCoordsFile(File f) {
coordsFile = f;
try {
@@ -823,7 +1022,7 @@ public void loadCoordsFile(File f) {
panel.setFont(font);
panel.setMoveStepSize(moveStepSize);
panel.setZoomStepSize(zoomStepSize);
- formatter = new FormatVertex(vertices, edgesio.getStats(),
+ formatter = new FormatVertex(vertices, edgesio.getLabels(),labelScale,0,0,0,0,false,edgesio.getStats(),
windowSizes, threads);
panel.setFormatter(formatter);
panel.fitData();
@@ -854,22 +1053,70 @@ public void loadCoordsFile(File f) {
*/
private JPanel buttonsPanel;
+ private void resetZooms()
+ {
+ panel.prepZoomPoint(false);
+ panel.prepZoomRegion(false);
+ zoomHighlighted.setSelected(false);
+ zoomRegion.setSelected(false);
+ }
+
+ private void inactiveButton(JToggleButton b)
+ {
+ //b.setContentAreaFilled(true);
+ b.setForeground(new JButton().getForeground());
+ b.setSelected(false);
+
+ }
+
+ private void setMode(modes mode)
+ {
+ if (mode==modes.handmode)
+ {
+ panel.handMode();
+ resetZooms();
+ magnifier.setSelected(false);
+ inactiveButton(magnifier);
+ //hand.setContentAreaFilled(false);
+ //hand.setOpaque(true);
+ hand.setForeground(Color.green);
+
+ }
+
+
+
+ if (mode==modes.magnifiermode)
+ {
+ panel.magnifierMode();
+ resetZooms();
+ hand.setSelected(false);
+ inactiveButton(hand);
+ //magnifier.setContentAreaFilled(false);
+ //magnifier.setOpaque(true);
+ magnifier.setForeground(Color.green);
+ }
+
+ if (mode==modes.nomode)
+ {
+ inactiveButton(hand);
+ inactiveButton(magnifier);
+ panel.noMode();
+
+ }
+
+
+ }
+
+
+
private void setButtons() {
buttonsPanel = new JPanel();
((FlowLayout) buttonsPanel.getLayout()).setAlignment(FlowLayout.LEFT);
container.add(buttonsPanel, BorderLayout.NORTH);
- JButton undo = new JButton("Undo");
- undo.addActionListener(new ActionListener() {
- public void actionPerformed(ActionEvent e) {
- panel.undo();
- panel.setPaintImage(); // TODO: move inside panel.xxxx()?
- panel.repaint();
- }
- });
- buttonsPanel.add(undo);
+
- JButton zoomIn = new JButton("In");
+ /*JButton zoomIn = new JButton("In");
zoomIn.addActionListener(new ActionListener() {
public void actionPerformed(ActionEvent e) {
VertexFitter f = new VertexFitter();
@@ -891,9 +1138,61 @@ public void actionPerformed(ActionEvent e) {
panel.repaint();
}
});
- buttonsPanel.add(zoomOut);
+ buttonsPanel.add(zoomOut);*/
+
+
+
+ //hand.setContentAreaFilled(false);
+ //hand.setOpaque(true);
+ //hand.setBackground(Color.black);
+ hand.setFocusable(false);
+
+ //hand.setContentAreaFilled(false);
+ hand.addActionListener(new ActionListener() {
+ public void actionPerformed(ActionEvent e) {
+ System.out.println(hand.getModel().isSelected());
+ if (hand.getModel().isSelected())
+ {
+ setMode(modes.handmode);
+ }
+ else
+ {
+ setMode(modes.nomode);
+
+ //hand.setSelected(false);
+ }
+
+ }
+ });
+ buttonsPanel.add(hand);
+
+
+
+ //magnifier.setContentAreaFilled(false);
+ //magnifier.setOpaque(true);
+ //magnifier.setBackground(Color.black);
+ magnifier.setFocusable(false);
+
+
+ //magnifier.setContentAreaFilled(false);
+ magnifier.addActionListener(new ActionListener() {
+ public void actionPerformed(ActionEvent e) {
+ System.out.println(magnifier.getModel().isSelected());
+ if (magnifier.getModel().isSelected())
+ {
+ setMode(modes.magnifiermode);
+ }
+ else
+ {
+ setMode(modes.nomode);
+ //hand.setSelected(false);
+ }
+
+ }
+ });
+ buttonsPanel.add(magnifier);
- JButton moveup = new JButton("Up");
+ /*JButton moveup = new JButton("Up");
moveup.addActionListener(new ActionListener() {
public void actionPerformed(ActionEvent e) {
panel.moveUp();
@@ -931,14 +1230,14 @@ public void actionPerformed(ActionEvent e) {
panel.repaint();
}
});
- buttonsPanel.add(moveright);
+ buttonsPanel.add(moveright);*/
JButton fit = new JButton("Fit");
fit.addActionListener(new ActionListener() {
public void actionPerformed(ActionEvent e) {
// TODO: SESS - twice?
panel.fitData();
- panel.fitData();
+ //panel.fitData();
System.out.println("Fit data still buggy.");
panel.setPaintImage(); // TODO: move inside panel.xxxx()?
panel.repaint();
@@ -996,8 +1295,35 @@ public void actionPerformed(ActionEvent e) {
panel.showVertices(false);
showVertices.setSelected(false);
buttonsPanel.add(showVertices);
+
+ JButton undo = new JButton("Undo");
+ undo.addActionListener(new ActionListener() {
+ public void actionPerformed(ActionEvent e) {
+ panel.undo();
+ panel.setPaintImage(); // TODO: move inside panel.xxxx()?
+ panel.repaint();
+ }
+ });
+ buttonsPanel.add(undo);
}
+ public void keyTyped(KeyEvent e) {
+ // Sys(e, "KEY TYPED: ");
+ }
+
+ /** Handle the key-pressed event from the text field. */
+ public void keyPressed(KeyEvent e) {
+ System.out.println(e.getKeyChar());
+ if (e.getKeyChar()==' ')
+ System.out.println("yes");
+ // displayInfo(e, "KEY PRESSED: ");
+ }
+
+ /** Handle the key-released event from the text field. */
+ public void keyReleased(KeyEvent e) {
+ // displayInfo(e, "KEY RELEASED: ");
+ }
+
class StyleHandler implements ItemListener {
public void itemStateChanged(ItemEvent e) {
int style = 0;
@@ -1037,6 +1363,15 @@ private boolean vertexCheck() {
return true;
}
+ public void componentHidden(ComponentEvent e) {}
+ public void componentMoved(ComponentEvent e) {}
+ public void componentShown(ComponentEvent e) {}
+
+ public void componentResized(ComponentEvent e) {
+ Dimension newSize = e.getComponent().getBounds().getSize();
+ System.out.println("size:"+newSize.width+" "+newSize.height);
+ }
+
}
// ////////////////////////////////////////////////////////////////
@@ -1095,6 +1430,8 @@ public String getExtension(File f) {
return ext;
}
+
+
}
// ////////////////////////////////////////////////////////////////
diff --git a/Java/src/Viewer2D/EdgesPanel.java b/Java/src/Viewer2D/EdgesPanel.java
index 404c3c20..a82bce42 100755
--- a/Java/src/Viewer2D/EdgesPanel.java
+++ b/Java/src/Viewer2D/EdgesPanel.java
@@ -20,20 +20,31 @@
package Viewer2D;
import java.awt.Color;
+import java.awt.Cursor;
import java.awt.Dimension;
import java.awt.Font;
+import java.awt.FontMetrics;
import java.awt.Graphics;
import java.awt.Graphics2D;
import java.awt.Rectangle;
import java.awt.RenderingHints;
+import java.awt.event.ActionListener;
+import java.awt.event.ActionEvent;
import java.awt.event.MouseEvent;
import java.awt.event.MouseListener;
import java.awt.event.MouseMotionListener;
+import java.awt.event.MouseWheelListener;
+
+import java.awt.event.MouseWheelEvent;
+
import java.awt.geom.Line2D;
import java.awt.geom.Rectangle2D;
+import java.awt.geom.Ellipse2D;
import java.awt.image.BufferedImage;
import java.io.BufferedWriter;
import java.io.File;
+import java.io.FileInputStream;
+import java.io.FileNotFoundException;
import java.io.FileWriter;
import java.io.IOException;
import java.text.MessageFormat;
@@ -42,12 +53,21 @@
import java.util.Map;
import javax.imageio.ImageIO;
+import javax.swing.JMenuItem;
import javax.swing.JOptionPane;
import javax.swing.JPanel;
+import javax.swing.JPopupMenu;
import javax.swing.JTextField;
import javax.swing.SwingUtilities;
import Jama.Matrix;
+import java.awt.BasicStroke;
+import Viewer2D.Label;
+import java.awt.AlphaComposite;
+import java.awt.FontFormatException;
+import java.awt.Point;
+import java.awt.Toolkit;
+import java.awt.Image;
/**
* SESS - 2014.05.13:
@@ -58,12 +78,14 @@
*
*/
public class EdgesPanel extends JPanel implements MouseListener,
- MouseMotionListener {
+ MouseMotionListener, MouseWheelListener {
private static final long serialVersionUID = -765273216020721560L;
private Edge[] edges;
private Vertex[] vertices;
+ private HashMap