001: import java.awt.*;
002: import java.awt.event.*;
003: import JSci.awt.*;
004: import JSci.maths.*;
005: import JSci.maths.matrices.*;
006: import JSci.maths.vectors.*;
007:
008: /**
009: * Wavefunction numerically solves the time-independent
010: * Schrödinger equation.
011: * @author Mark Hale
012: * @version 1.2
013: */
014: public final class Wavefunction extends Frame implements Runnable {
015: private final Runnable animator = this ;
016: private final int N = 200;
017: private final LineGraph graph;
018: private final GraphModel model = new GraphModel();
019: private final Label kLabel = new Label("Harmonic coupling");
020: private final TextField kField = new TextField("0.0", 5);
021: private final Label hLabel = new Label("Anharmonic coupling");
022: private final TextField hField = new TextField("0.0", 5);
023: private final Button probButton = new Button("Probability");
024: private final Button animButton = new Button("Evolve");
025: private final Label statusLabel = new Label("-", Label.CENTER);
026: private final Button incButton = new Button("+");
027: private final Button decButton = new Button("-");
028: private volatile Thread animateThread = null;
029: /**
030: * Eigenstates.
031: */
032: private final AbstractDoubleVector eigenstates[] = new AbstractDoubleVector[N];
033: /**
034: * Eigenvalues.
035: */
036: private double eigenvalues[];
037: /**
038: * Harmonic coupling.
039: */
040: private double k = 0.0;
041: /**
042: * Anharmonic coupling.
043: */
044: private double h = 0.0;
045: /**
046: * Kinetic energy operator.
047: */
048: private final DoubleTridiagonalMatrix T = new DoubleTridiagonalMatrix(
049: N);
050: /**
051: * Potential energy operator.
052: */
053: private DoubleDiagonalMatrix V;
054:
055: /**
056: * Runs an instance of Wavefunction.
057: */
058: public static void main(String arg[]) {
059: new Wavefunction();
060: }
061:
062: private static void setDefaultSize(Component c, int defaultWidth,
063: int defaultHeight) {
064: Dimension screenSize = Toolkit.getDefaultToolkit()
065: .getScreenSize();
066: final int width = (defaultWidth < screenSize.width) ? defaultWidth
067: : screenSize.width;
068: final int height = (defaultHeight < screenSize.height) ? defaultHeight
069: : screenSize.height;
070: c.setSize(width, height);
071: }
072:
073: /**
074: * Constructs Wavefunction.
075: */
076: public Wavefunction() {
077: super ("Wavefunction");
078: constructKineticEnergyOperator();
079: calculateHamiltonian();
080: graph = new LineGraph(model);
081: graph.setYExtrema(-0.5f, 0.5f);
082: // register listeners
083: addWindowListener(new WindowAdapter() {
084: public void windowClosing(WindowEvent evt) {
085: dispose();
086: System.exit(0);
087: }
088: });
089: incButton.addActionListener(new ActionListener() {
090: public void actionPerformed(ActionEvent evt) {
091: model.incrementLevel();
092: updateStatusLabel();
093: }
094: });
095: decButton.addActionListener(new ActionListener() {
096: public void actionPerformed(ActionEvent evt) {
097: model.decrementLevel();
098: updateStatusLabel();
099: }
100: });
101: kField.addActionListener(new ActionListener() {
102: public void actionPerformed(ActionEvent evt) {
103: k = Double.valueOf(kField.getText()).doubleValue();
104: calculateHamiltonian();
105: updateStatusLabel();
106: }
107: });
108: hField.addActionListener(new ActionListener() {
109: public void actionPerformed(ActionEvent evt) {
110: h = Double.valueOf(hField.getText()).doubleValue();
111: calculateHamiltonian();
112: updateStatusLabel();
113: }
114: });
115: probButton.addActionListener(new ActionListener() {
116: public void actionPerformed(ActionEvent evt) {
117: if (model.isShowingProbability()) {
118: model.showAmplitude();
119: graph.setYExtrema(-0.5f, 0.5f);
120: probButton.setLabel("Probability");
121: } else {
122: model.showProbability();
123: graph.setYExtrema(0.0f, 0.25f);
124: probButton.setLabel("Amplitude");
125: }
126: }
127: });
128: animButton.addActionListener(new ActionListener() {
129: public void actionPerformed(ActionEvent evt) {
130: if (animateThread == null) {
131: animateThread = new Thread(animator);
132: animateThread.start();
133: animButton.setLabel("Freeze");
134: } else {
135: animateThread = null;
136: animButton.setLabel("Evolve");
137: }
138: }
139: });
140: // construct gui
141: Panel levelPanel = new Panel();
142: levelPanel.setLayout(new GridLayout(2, 1));
143: levelPanel.add(incButton);
144: levelPanel.add(decButton);
145: Panel optPanel = new Panel();
146: optPanel.add(kLabel);
147: optPanel.add(kField);
148: optPanel.add(hLabel);
149: optPanel.add(hField);
150: optPanel.add(probButton);
151: optPanel.add(animButton);
152: add(statusLabel, "North");
153: add(graph, "Center");
154: add(levelPanel, "East");
155: add(optPanel, "South");
156: setDefaultSize(this , 650, 450);
157: updateStatusLabel();
158: setVisible(true);
159: }
160:
161: public void run() {
162: while (animateThread == Thread.currentThread()) {
163: EventQueue.invokeLater(new Runnable() {
164: public void run() {
165: model.evolve();
166: }
167: });
168: try {
169: Thread.sleep(100);
170: } catch (InterruptedException e) {
171: }
172: }
173: }
174:
175: /**
176: * Updates the status label.
177: */
178: private void updateStatusLabel() {
179: final int level = model.getLevel();
180: statusLabel.setText("Energy [" + level + "] = "
181: + (float) eigenvalues[level]);
182: }
183:
184: /**
185: * A custom graph model class.
186: * This provides a graph interface onto the eigenstates
187: * and eigenvalues.
188: * There are various methods to alter the data view.
189: */
190: private final class GraphModel extends AbstractGraphModel implements
191: Graph2DModel {
192: /**
193: * Time.
194: */
195: private double t = 0.0;
196: /**
197: * Energy level.
198: */
199: private int level = 0;
200: /**
201: * Show probability or amplitude.
202: */
203: private boolean showProb = false;
204: /**
205: * Series.
206: */
207: private final static int SERIES_WAVEFUNCTION = 0;
208: private final static int SERIES_POTENTIAL = 1;
209: private int series = SERIES_WAVEFUNCTION;
210:
211: public GraphModel() {
212: }
213:
214: public float getXCoord(int i) {
215: return i * 2.0f / (N - 1) - 1.0f;
216: }
217:
218: public float getYCoord(int i) {
219: if (series == SERIES_WAVEFUNCTION) {
220: final double amp = eigenstates[level].getComponent(i);
221: if (showProb)
222: return (float) (amp * amp);
223: else
224: return (float) (amp * Math.cos(eigenvalues[level]
225: * t));
226: } else if (series == SERIES_POTENTIAL)
227: return (float) (V.getElement(i, i) - eigenvalues[level]);
228: else
229: return 0.0f;
230: }
231:
232: public void resetTime() {
233: t = 0.0;
234: fireGraphSeriesUpdated(SERIES_WAVEFUNCTION);
235: }
236:
237: public void evolve() {
238: t += 2.0;
239: fireGraphSeriesUpdated(SERIES_WAVEFUNCTION);
240: }
241:
242: public void incrementLevel() {
243: if (level < N - 1)
244: level++;
245: fireGraphSeriesUpdated(series);
246: }
247:
248: public void decrementLevel() {
249: if (level > 0)
250: level--;
251: fireGraphSeriesUpdated(series);
252: }
253:
254: public int getLevel() {
255: return level;
256: }
257:
258: public void showAmplitude() {
259: showProb = false;
260: fireGraphSeriesUpdated(SERIES_WAVEFUNCTION);
261: }
262:
263: public void showProbability() {
264: showProb = true;
265: fireGraphSeriesUpdated(SERIES_WAVEFUNCTION);
266: }
267:
268: public boolean isShowingAmplitude() {
269: return !showProb;
270: }
271:
272: public boolean isShowingProbability() {
273: return showProb;
274: }
275:
276: public int seriesLength() {
277: return N;
278: }
279:
280: public void firstSeries() {
281: series = SERIES_WAVEFUNCTION;
282: }
283:
284: public boolean nextSeries() {
285: series++;
286: if (series > 1)
287: return false;
288: else
289: return true;
290: }
291: }
292:
293: // Solve the Schrodinger equation
294:
295: /**
296: * Constructs the kinetic energy operator.
297: */
298: private void constructKineticEnergyOperator() {
299: T.setElement(0, 0, 2.0);
300: T.setElement(0, 1, -1.0);
301: int i = 1;
302: for (; i < N - 1; i++) {
303: T.setElement(i, i - 1, -1.0);
304: T.setElement(i, i, 2.0);
305: T.setElement(i, i + 1, -1.0);
306: }
307: T.setElement(i, i - 1, -1.0);
308: T.setElement(i, i, 2.0);
309: }
310:
311: /**
312: * Constructs the Hamiltonian
313: * and solves the Schrödinger equation.
314: */
315: private void calculateHamiltonian() {
316: // potential energy operator
317: final double pot[] = new double[N];
318: for (int i = 0; i < N; i++)
319: pot[i] = potential(i);
320: V = new DoubleDiagonalMatrix(pot);
321: // Hamiltonian
322: final AbstractDoubleSquareMatrix H = T.add(V);
323: // solve
324: try {
325: eigenvalues = LinearMath
326: .eigenSolveSymmetric(H, eigenstates);
327: } catch (MaximumIterationsExceededException e) {
328: System.err.println(e.getMessage());
329: System.exit(-1);
330: }
331: sortStates();
332: model.resetTime();
333: }
334:
335: /**
336: * Potential energy function.
337: */
338: private double potential(int i) {
339: final double x = model.getXCoord(i);
340: return k * x * x + h * x * x * x * x;
341: }
342:
343: /**
344: * Sorts the eigenstates and eigenvalues.
345: * (Bidirectional bubble sort.)
346: */
347: private void sortStates() {
348: int i;
349: int limit = eigenstates.length;
350: int st = -1;
351: AbstractDoubleVector vec;
352: double val;
353: while (st < limit) {
354: boolean flipped = false;
355: st++;
356: limit--;
357: for (i = st; i < limit; i++) {
358: if (eigenvalues[i] > eigenvalues[i + 1]) {
359: vec = eigenstates[i];
360: val = eigenvalues[i];
361: eigenstates[i] = eigenstates[i + 1];
362: eigenvalues[i] = eigenvalues[i + 1];
363: eigenstates[i + 1] = vec;
364: eigenvalues[i + 1] = val;
365: flipped = true;
366: }
367: }
368: if (!flipped)
369: return;
370: for (i = limit; --i >= st;) {
371: if (eigenvalues[i] > eigenvalues[i + 1]) {
372: vec = eigenstates[i];
373: val = eigenvalues[i];
374: eigenstates[i] = eigenstates[i + 1];
375: eigenvalues[i] = eigenvalues[i + 1];
376: eigenstates[i + 1] = vec;
377: eigenvalues[i + 1] = val;
378: flipped = true;
379: }
380: }
381: if (!flipped)
382: return;
383: }
384: }
385: }
|