In [1]:
import matplotlib.pyplot as plt
import os
import pandas as pd
import numpy as np
import quantumpropagator as qp
from quantumpropagator import readWholeH5toDict, abs2
#from __future__ import print_function
from ipywidgets import interact#, interactive, fixed, interact_manual
import ipywidgets as widgets

def heatThis(h5File,gams,thes,vmaxV,save=None):
    save = save or False
    from mpl_toolkits.axes_grid1 import make_axes_locatable
    wf   = qp.retrieve_hdf5_data(h5File,'WF')
    time = qp.retrieve_hdf5_data(h5File,'Time')[0]
    fig = plt.figure(figsize=(10, 10), dpi= 80, facecolor='w', edgecolor='k')
    plt.title('Time = {:10.5f} fs'.format(time))
    plt.xlabel('Gamma')
    plt.ylabel('Theta')
    
    # this is to get a nice colorbar on the side
    ax = plt.gca()
    aaa = np.rad2deg
    ext = [aaa(gams[0]),aaa(gams[-1]),aaa(thes[0])*2,aaa(thes[-1])*2]
    im = ax.imshow(qp.abs2(wf), extent=ext, cmap='hot', vmax=vmaxV)
    #im = ax.imshow(qp.abs2(wf), cmap='PuBu_r', vmax=0.4)
    divider = make_axes_locatable(ax)
    cax = divider.append_axes("right", size="5%", pad=0.05)

    plt.colorbar(im, cax=cax)
    
    if save:
        fig.savefig(os.path.splitext(h5File)[0] + '.png')
        plt.close()
    
def f(x):
    return heatThis(os.path.join(a,filesList[x]),gams,thes,vmaxV)

In [2]:
subfolders = sorted([ dir for dir in os.listdir('.') if os.path.isdir(dir) and dir != '.ipynb_checkpoints'])
print(''.join(['{} -> {}\n'.format(a,b) for a,b in enumerate(subfolders)]))


0 -> input_FinerGrid_0000
1 -> input_FinerGrid_0001
2 -> input_FinerGrid_0002
3 -> input_FinerGrid_0003
4 -> input_FinerGrid_0004
5 -> input_FinerGrid_0005


In [3]:
a = subfolders[-1]

In [4]:
filesList = [ fn for fn in sorted(os.listdir(a)) if fn[:8] == 'Gaussian' and fn[-3:] == '.h5']
outh5 = os.path.join(a,'allInput.h5')
dictio = readWholeH5toDict(outh5)
gams,thes = dictio['gams'],dictio['thes']
lastV = len(filesList)-1
# dictio.keys()
zeroWF = qp.retrieve_hdf5_data(os.path.join(a,filesList[0]),'WF')
vmaxV = abs2(zeroWF).max()

In [5]:
interact(f, x = widgets.IntSlider(min=0,max=lastV,step=1,value=lastV));
createimages = False

if createimages:
    for filn in filesList:
        filna = os.path.join(a,filn)
        heatThis(filna, gams, thes, vmaxV, save=True)


1d slices


In [6]:
gamL,theL = (qp.retrieve_hdf5_data(os.path.join(a,filesList[0]),'WF')).shape
gamsT = np.rad2deg(gams)
thesT = np.rad2deg(thes)*2

In [7]:
def sliceGammas(h5File, gamma, gamsT):
    wf  = qp.retrieve_hdf5_data(h5File,'WF')
    fig = plt.figure(figsize=(11, 6), dpi= 80, facecolor='w', edgecolor='k')
    ys = wf[gamma,:]
    time = qp.retrieve_hdf5_data(h5File,'Time')[0]
    plt.title('Time = {:10.5f} fs --- gammaL = {:8.3f}'.format(time,gamsT[gamma]))
    plt.ylim(-0.5,0.5)
    plt.plot(thesT, np.real(ys), linewidth=1,ls='--')
    plt.plot(thesT, np.imag(ys), linewidth=1,ls='--')
    plt.plot(thesT, abs2(ys), linewidth=3,ls='-')

def sliceThetas(h5File, theta):
    wf  = qp.retrieve_hdf5_data(h5File,'WF')
    fig = plt.figure(figsize=(11, 6), dpi= 80, facecolor='w', edgecolor='k')
    ys = wf[:,theta]
    time = qp.retrieve_hdf5_data(h5File,'Time')[0]
    plt.title('Time = {:10.5f} fs --- thetaL = {:8.3f}'.format(time,thesT[theta]))
    plt.ylim(-0.5,0.5)
    plt.plot(gamsT, np.real(ys), linewidth=1,ls='--')
    plt.plot(gamsT, np.imag(ys), linewidth=1,ls='--')
    plt.plot(gamsT, abs2(ys), linewidth=3,ls='-')    
    
def fgam(file_number,gamma):
    return sliceGammas(os.path.join(a, filesList[file_number]), gamma, gamsT)

def fthe(file_number,theta):
    return sliceThetas(os.path.join(a, filesList[file_number]), theta)

In [8]:
interact(fgam, file_number = widgets.IntSlider(min=0,max=lastV,step=1,value=lastV), gamma = widgets.IntSlider(min=0,max=gamL-1,step=1,value=10));



In [9]:
interact(fthe, file_number = widgets.IntSlider(min=0,max=lastV,step=1,value=lastV), theta = widgets.IntSlider(min=0,max=theL-1,step=1,value=16));



In [10]:
outfn = os.path.join(a,'output')

In [11]:
data = pd.read_csv(outfn, delim_whitespace=True, header=None);
data.columns = ['steps','fs','Norm Deviation','Kinetic','Potential','Total','Total deviation']
data


Out[11]:
steps fs Norm Deviation Kinetic Potential Total Total deviation
0 0 0.000000 0.396577 0.000860 0.001375 0.002235 0.000000e+00
1 100 1.000993 0.396581 0.000968 0.001268 0.002235 6.454296e-08
2 200 2.001985 0.396594 0.001185 0.001050 0.002235 2.394490e-07
3 300 3.002978 0.396612 0.001397 0.000838 0.002235 4.882732e-07
4 400 4.003971 0.396634 0.001519 0.000716 0.002235 7.740499e-07
5 500 5.004964 0.396657 0.001529 0.000705 0.002234 1.067840e-06
6 600 6.005956 0.396679 0.001441 0.000793 0.002234 1.348027e-06
7 700 7.006949 0.396700 0.001285 0.000949 0.002234 1.597554e-06
8 800 8.007942 0.396717 0.001102 0.001132 0.002234 1.802655e-06
9 900 9.008935 0.396730 0.000942 0.001291 0.002233 1.953766e-06
10 1000 10.009927 0.396737 0.000848 0.001386 0.002233 2.044166e-06
11 1100 11.010920 0.396740 0.000840 0.001393 0.002233 2.070173e-06
12 1200 12.011913 0.396737 0.000920 0.001313 0.002233 2.031568e-06
13 1300 13.012906 0.396729 0.001068 0.001165 0.002233 1.930433e-06
14 1400 14.013898 0.396716 0.001249 0.000985 0.002234 1.772077e-06
15 1500 15.014891 0.396700 0.001413 0.000821 0.002234 1.564396e-06
16 1600 16.015884 0.396680 0.001518 0.000716 0.002234 1.319128e-06
17 1700 17.016877 0.396659 0.001525 0.000709 0.002234 1.051455e-06
18 1800 18.017869 0.396638 0.001424 0.000810 0.002235 7.824501e-07
19 1900 19.018862 0.396619 0.001249 0.000986 0.002235 5.406982e-07
20 2000 20.019855 0.396605 0.001081 0.001154 0.002235 3.588816e-07
21 2100 21.020847 0.396597 0.000984 0.001251 0.002235 2.630683e-07
22 2200 22.021840 0.396595 0.000986 0.001250 0.002235 2.638712e-07
23 2300 23.022833 0.396601 0.001066 0.001169 0.002235 3.556312e-07
24 2400 24.023826 0.396613 0.001200 0.001035 0.002235 5.241483e-07
25 2500 25.024818 0.396631 0.001346 0.000888 0.002235 7.504349e-07
26 2600 26.025811 0.396651 0.001452 0.000783 0.002234 1.010322e-06
27 2700 27.026804 0.396672 0.001476 0.000758 0.002234 1.279160e-06
28 2800 28.027797 0.396693 0.001412 0.000822 0.002234 1.535738e-06
29 2900 29.028789 0.396713 0.001283 0.000951 0.002234 1.762559e-06
... ... ... ... ... ... ... ...
71 7100 71.070484 0.396667 0.001446 0.000788 0.002234 1.262035e-06
72 7200 72.071477 0.396688 0.001354 0.000880 0.002234 1.505815e-06
73 7300 73.072470 0.396707 0.001215 0.001019 0.002234 1.716856e-06
74 7400 74.073462 0.396722 0.001068 0.001165 0.002234 1.883971e-06
75 7500 75.074455 0.396734 0.000951 0.001282 0.002233 1.999392e-06
76 7600 76.075448 0.396740 0.000895 0.001338 0.002233 2.058166e-06
77 7700 77.076441 0.396742 0.000910 0.001323 0.002233 2.057938e-06
78 7800 78.077433 0.396738 0.000993 0.001241 0.002233 1.998805e-06
79 7900 79.078426 0.396729 0.001123 0.001110 0.002234 1.883522e-06
80 8000 80.079419 0.396716 0.001269 0.000964 0.002234 1.717429e-06
81 8100 81.080412 0.396699 0.001384 0.000849 0.002234 1.508950e-06
82 8200 82.081404 0.396679 0.001436 0.000798 0.002234 1.271306e-06
83 8300 83.082397 0.396658 0.001415 0.000820 0.002234 1.022659e-06
84 8400 84.083390 0.396639 0.001334 0.000901 0.002235 7.845924e-07
85 8500 85.084383 0.396621 0.001216 0.001019 0.002235 5.773707e-07
86 8600 86.085375 0.396608 0.001088 0.001147 0.002235 4.218224e-07
87 8700 87.086368 0.396600 0.001002 0.001233 0.002235 3.380606e-07
88 8800 88.087361 0.396599 0.001007 0.001228 0.002235 3.399480e-07
89 8900 89.088354 0.396605 0.001104 0.001131 0.002235 4.281090e-07
90 9000 90.089346 0.396617 0.001250 0.000985 0.002235 5.902243e-07
91 9100 91.090339 0.396633 0.001391 0.000844 0.002235 8.053670e-07
92 9200 92.091332 0.396653 0.001474 0.000760 0.002234 1.048588e-06
93 9300 93.092324 0.396674 0.001477 0.000757 0.002234 1.296659e-06
94 9400 94.093317 0.396694 0.001395 0.000839 0.002234 1.530236e-06
95 9500 95.094310 0.396712 0.001256 0.000978 0.002234 1.733677e-06
96 9600 96.095303 0.396728 0.001099 0.001134 0.002234 1.895023e-06
97 9700 97.096295 0.396738 0.000961 0.001273 0.002233 2.005354e-06
98 9800 98.097288 0.396744 0.000879 0.001355 0.002233 2.059684e-06
99 9900 99.098281 0.396745 0.000878 0.001355 0.002233 2.055275e-06
100 9999 100.089264 0.396741 0.000956 0.001278 0.002233 1.992293e-06

101 rows × 7 columns


In [12]:
data.plot(title = 'Norm Deviation', x='fs', y = 'Norm Deviation', figsize=(15,4));



In [13]:
data.plot(title = 'Kinetic', x='fs', y = 'Kinetic', figsize=(15,4));



In [14]:
data.plot(title = 'Potential', x='fs', y = 'Potential', figsize=(15,4));



In [15]:
data['Kinetic_Moved'] = data['Kinetic'] + data['Potential'][0]
data.plot(title = 'Comparison Potential Total Kinetic', x=['fs'] ,y=['Kinetic_Moved','Potential','Total'], figsize=(15,5));


Expected values


In [9]:
def expected(h5file):
    wf  = qp.retrieve_hdf5_data(h5file,'WF')
    wf = wf /np.linalg.norm(wf)
    wfA = abs2(wf)
    gamAvg = np.sum(wfA,axis=1)
    theAvg = np.sum(wfA,axis=0)
    gamA = sum([ gamsT[i] * gamAvg[i] for i in range(gamAvg.size) ])
    theA = sum([ thesT[i] * theAvg[i] for i in range(theAvg.size) ])
    print('Expected Values:\nTheta: {:10.3f} \nGamma: {:10.3f}'.format(theA, gamA))

def fexp(x):
    return expected(os.path.join(a,filesList[x]))

In [10]:
interact(fexp, x = widgets.IntSlider(min=0,max=lastV,step=1,value=lastV));