In [1]:
#------------------------------------------------------------------------------------
import numpy as np
import pycuda.gpuarray as gpuarray
from pycuda.tools import make_default_context
import matplotlib as matplotlib
import pylab as plt
from mpl_toolkits.mplot3d import Axes3D
#-------------------------------------------------------------------------------------
from pywignercuda_path import SetPyWignerCUDA_Path
SetPyWignerCUDA_Path()
from GPU_Dirac3D import *
In [2]:
#import pyqtgraph as pg
In [3]:
%matplotlib inline
In [4]:
class DefinePotential_Dirac3D( GPU_Dirac3D ):
def __init__ (self):
#....................Defining the geometry.....................................
gridDIM_X = 256
gridDIM_Y = 256
gridDIM_Z = 256
X_amplitude =15
Y_amplitude =15
Z_amplitude =15
dt = 0.01
timeSteps = 1200
skipFrames = 100
#...................Defining the kinematic-dynamical constants.................
self.mass = 1.
self.c = 1.
self.hBar = 1.
#...................Defining the potential and initial state parameters........
hBar = 1
self.px = 0.
self.py = 5.
self.pz = 0.
self.energy = self.c*np.sqrt( (self.mass*self.c)**2 + self.px**2 + self.py**2 + self.pz**2 )
#...................Defining the potential.....................................
V0 = 7.
w = 0.3
self.Potential_0_String = '{0}*0.5*( 1. + tanh((y-2.)/{1}) )'.format(V0,w)
self.Potential_1_String = '0.'
self.Potential_2_String = '0.'
self.Potential_3_String = '0.'
#..................Defining the output directory/file .........................
self.fileName = '/home/rcabrera/DATA/Dirac3D/Klein_Y.hdf5'
self.Compute_Ehrenfest_P = True
self.antiParticleStepFiltering = False
#..............................................................................
amplitude = (X_amplitude,Y_amplitude,Z_amplitude)
gridDIM = (gridDIM_X,gridDIM_Y,gridDIM_Z)
self.Compute_Ehrenfest_P = True
GPU_Dirac3D.__init__( self,gridDIM, amplitude ,dt,timeSteps, skipFrames = skipFrames,frameSaveMode='Spinor')
def Set_Initial_Condition_SpinUp (self):
def gaussian(x,y,z):
y0 = -3
return np.exp( - ( x**2 + (y-y0)**2 + z**2)/2. )
self.Psi_init = self.Spinor_Particle_SpinUp( (self.px,self.py,self.pz) , gaussian )
norm = self.Norm( self.Psi_init )
#self.FilterElectrons(1)
self.Psi_init /= norm
"""def Set_Initial_Condition_SpinDown (self):
def gaussian(x,y,z):
return np.exp( - ( x**2 + y**2 + z**2)/2. )
self.Psi_init = self.Spinor_Particle_SpinDown( (self.px,self.py,self.pz) , gaussian )
norm = self.Norm_X_CPU( self.Psi_init )
self.Psi_init /= norm """
In [5]:
instance = DefinePotential_Dirac3D()
In [6]:
instance = DefinePotential_Dirac3D()
instance.Set_Initial_Condition_SpinUp()
#instance_Dirac3D.FilterElectrons(1)
print ' '
print ' dx = ', instance.dX
Psi_end = instance.Run ()
In [7]:
instance.Psi_end.shape
Out[7]:
In [8]:
axis_font = {'size':'24'}
fig, ax = plt.subplots(figsize=(20, 7))
ax.plot( instance.timeRange , instance.X_average , 'g',label= '$x^1 $')
ax.plot( instance.timeRange , instance.Y_average , 'r',label= '$x^2 $')
ax.plot( instance.timeRange , instance.Z_average , 'b--',label= '$x^3 $')
ax.set_xlabel(r'$t$',**axis_font)
ax.legend(bbox_to_anchor=(1.05, 0.5), loc=1, prop={'size':22})
Out[8]:
In [9]:
if instance.Compute_Ehrenfest_P == True:
axis_font = {'size':'24'}
fig, ax = plt.subplots(figsize=(20, 7))
ax.plot( instance.timeRange , instance.Px_average , 'g',label= '$p^1 $')
ax.plot( instance.timeRange , instance.Py_average , 'r',label= '$p^2 $')
ax.plot( instance.timeRange , instance.Pz_average , 'b--',label= '$p^3 $')
ax.set_ylim(-6,6)
ax.set_xlabel(r'$t$',**axis_font)
ax.legend(bbox_to_anchor=(1.05, 0.5), loc=1, prop={'size':22})
In [10]:
axis_font = {'size':'24'}
fig, ax = plt.subplots(figsize=(20, 7))
ax.plot( instance.timeRange , np.gradient( instance.X_average , instance.dt) , 'g',
label= '$\\frac{dx^1}{dt} $')
ax.plot( instance.timeRange , instance.Alpha1_average ,'r--' ,label='$c \\alpha^1$')
ax.set_xlabel(r'$t$',**axis_font)
ax.legend(bbox_to_anchor=(1.05, 0.5), loc=1, prop={'size':22})
Out[10]:
In [11]:
axis_font = {'size':'24'}
fig, ax = plt.subplots(figsize=(20, 7))
ax.plot( instance.timeRange , np.gradient( instance.Y_average , instance.dt) , 'g',
label= '$\\frac{dy^1}{dt} $')
ax.plot( instance.timeRange , instance.Alpha2_average ,'r--' ,label='$c \\alpha^2$')
ax.set_xlabel(r'$t$',**axis_font)
ax.legend(bbox_to_anchor=(1.05, 0.5), loc=1, prop={'size':22})
Out[11]:
In [12]:
axis_font = {'size':'24'}
fig, ax = plt.subplots(figsize=(20, 7))
ax.plot( instance.timeRange , np.gradient( instance.Z_average , instance.dt) , 'g',
label= '$\\frac{dx^3}{dt} $')
ax.plot( instance.timeRange , instance.Alpha3_average ,'r--' ,label='$c \\alpha^3$')
ax.set_ylim(-0.1,0.1)
ax.set_xlabel(r'$t$',**axis_font)
ax.legend(bbox_to_anchor=(1.05, 0.5), loc=1, prop={'size':22})
Out[12]:
In [13]:
if instance.Compute_Ehrenfest_P == True:
axis_font = {'size':'24'}
fig, ax = plt.subplots(figsize=(20, 7))
ax.plot( instance.timeRange , np.gradient( instance.Px_average , instance.dt) , 'g',
label= '$\\frac{dp^1}{dt} $')
ax.set_xlabel(r'$t$',**axis_font)
ax.legend(bbox_to_anchor=(1.05, 0.5), loc=1, prop={'size':22})
In [14]:
if instance.Compute_Ehrenfest_P == True:
axis_font = {'size':'24'}
fig, ax = plt.subplots(figsize=(20, 7))
ax.plot( instance.timeRange , np.gradient( instance.Py_average , instance.dt) , 'g',
label= '$\\frac{dp^2}{dt} $')
ax.set_xlabel(r'$t$',**axis_font)
ax.legend(bbox_to_anchor=(1.05, 0.5), loc=1, prop={'size':22})
In [15]:
if instance.Compute_Ehrenfest_P == True:
axis_font = {'size':'24'}
fig, ax = plt.subplots(figsize=(20, 7))
ax.plot( instance.timeRange , np.gradient( instance.Pz_average , instance.dt) , 'g',
label= '$\\frac{dp^3}{dt} $')
ax.set_xlabel(r'$t$',**axis_font)
ax.legend(bbox_to_anchor=(1.05, 0.5), loc=1, prop={'size':22})
In [16]:
if instance.Compute_Ehrenfest_P == True:
axis_font = {'size':'24'}
fig, ax = plt.subplots(figsize=(20, 7))
ax.plot( instance.timeRange , instance.K_Energy_average + instance.Potential_0_average , 'g',
label= '$ Energy$')
ax.set_ylim(0,7)
ax.set_xlabel(r'$t$',**axis_font)
ax.grid('on')
ax.legend(bbox_to_anchor=(1.05, 0.5), loc=1, prop={'size':22})
In [22]:
fig, ax = plt.subplots(figsize=(20, 10))
#ax = fig.gca(projection='3d')
ax.scatter(
instance.X_average,
instance.Y_average )
ax.set_xlabel(r'$x$',**axis_font)
ax.set_ylabel(r'$y$',**axis_font)
ax.set_xlim(-1. ,1.)
ax.set_ylim(-5 ,5)
#ax.set_zlim(-5,0)
Out[22]:
In [23]:
fig, ax = plt.subplots(figsize=(20, 10))
ax = fig.gca(projection='3d')
ax.scatter(
instance.X_average,
instance.Y_average,
instance.Z_average
)
ax.set_xlabel(r'$x$',**axis_font)
ax.set_ylabel(r'$y$',**axis_font)
ax.set_zlabel(r'$z$',**axis_font)
ax.set_ylim(-2 ,6)
ax.set_zlim(-0.1,0.1)
ax.set_xlim(-1,1)
Out[23]:
In [ ]: