diff --git a/QuantumMovement/crankNicolson/animate.py b/QuantumMovement/crankNicolson/animate.py index 943cc959..cf33a8d0 100644 --- a/QuantumMovement/crankNicolson/animate.py +++ b/QuantumMovement/crankNicolson/animate.py @@ -1,7 +1,15 @@ +# ! ! !TO DO +# Fix sacle and cutoff to work with new implementation of changing between mod2/phase/real/imag representation +# Ek: Vermell, V: Blau, E: Lila import matplotlib from matplotlib import animation import matplotlib.pyplot as plt import matplotlib.patheffects as peffects +from matplotlib.colors import LinearSegmentedColormap +#basic_cols=['#75b765', '#101010', '#ffd700'] +basic_cols=['#00ffff', '#101010', '#ff0000'] +div_cmap=LinearSegmentedColormap.from_list('bl_div_cmap', basic_cols) + if (__name__ == "crankNicolson.animate"): import crankNicolson.crankNicolson2D as mathPhysics from kivy.garden.matplotlib.backend_kivyagg import FigureCanvasKivyAgg, NavigationToolbar2Kivy @@ -75,6 +83,57 @@ def drawOptimized(self): # fig.canvas.mpl_disconnect(fig.canvas.manager.key_press_handler_id) +def cleanFigClose(fig): + """ + Close figure trying to avoid memory leaks!!! + :param fig: figure to be closed cleanly + """ + try: fig.canvas.canvas.clear() # this is for the kivy widget. Makes no sense otherwise + except: pass + fig.clf() + plt.close(fig) + +#from colorsys import hls_to_rgb +from matplotlib.colors import hsv_to_rgb +def Complex2HSV(z, rmin=None, rmax=None, hue_start=90): + """ + https://stackoverflow.com/a/36082859 + """ + # get amplidude of z and limit to [rmin, rmax] + amp = np.abs(z) + #amp = mathPhysics.abs2(z) + if rmin != None: amp = np.where(amp < rmin, rmin, amp) + if rmax != None: amp = np.where(amp > rmax, rmax, amp) + ph = np.angle(z, deg=True) + hue_start + # HSV are values in range [0,1] + h = (ph % 360) / 360 + s = 0.85 * np.ones_like(h) + v = amp/np.max(amp)#(amp -rmin) / (rmax - rmin) + return hsv_to_rgb(np.dstack((h,s,v))) + +# We transpose!!! Because we take [i,j] -> i corresponds to x, j to y +# Whereas in matrix representation (imshow), i would be row (y) and j column (x) +def plotComplexData(arr, type='mod2'): + if type=='phase': + return Complex2HSV(arr.T) + elif type=='real': + return arr.real.T + elif type=='imag': + return arr.imag.T + return mathPhysics.abs2(arr).T + +def plotComplex(arr, ax, type, range): + if type=='phase': + return ax.imshow(plotComplexData(arr, type), origin='lower', extent=range, aspect='equal', interpolation='none') + elif type=='mod2': + return ax.imshow(plotComplexData(arr, type), origin='lower', extent=range, aspect='equal', + cmap="viridis", interpolation='none') + else: + maxV = np.max(plotComplexData(arr, type)) + minV = np.min(plotComplexData(arr, type)) + sym = max(maxV, abs(minV)) + return ax.imshow(plotComplexData(arr, type), origin='lower', extent=range, aspect='equal', + cmap=div_cmap, interpolation='none', vmin=-sym,vmax=sym) optimizeGraphics = True @@ -82,12 +141,13 @@ class QuantumAnimation: # inches def __init__(self, QSystem, width=6.4, height=4.8, dtSim=0.01, dtAnim=0.05, duration=None, realTime=True, showEnergy=False, forceEnergy=False, showNorm=False, forceNorm=False, showMomentum=False, showPotential=False, varyingPotential=False, + psiRepresentation="mod2", momRepresentation="mod2", potentialCutoff=None, psiCutoff=None, scalePsi=False, scaleMom=False, zoomMom=1., scalePot=True, extraCommands=[], extraUpdates=None, extraUpdatesStart=False, extraUpdatesUpdate=False, isKivy=False, stepsPerFrame=0, debugTime=False, callbackProgress=False, isFocusable=True, drawClassical=False, drawClassicalTrace=False, drawExpected=False, drawExpectedTrace=False, unit_dist='Å', unit_time='fs', unit_energy='eV', unit_mom=r'$\frac{1}{2}\hbar Å^{-1}$', - toolbar=False, customPlot=None, customPlotUpdate=False, customPlotFull=False): + toolbar=False, customPlot=None, customPlotUpdate=False, customPlotFull=False,): """ Declaration of a Quantum Animation :param QSystem: Quantum System. The physics happen here, see crankNikolson2D QuantumSystem class. @@ -144,6 +204,9 @@ def __init__(self, QSystem, width=6.4, height=4.8, self.showPotential = showPotential self.varyingPotential = varyingPotential self.potentialCutoff = potentialCutoff + + self.psiRepresentation = psiRepresentation + self.momRepresentation = momRepresentation self.psiCutoff = psiCutoff self.scalePsi = scalePsi self.scaleMom = scaleMom @@ -247,7 +310,7 @@ def __init__(self, QSystem, width=6.4, height=4.8, def reset_plot(self, width=None, height=None, showEnergy=None, forceEnergy=None, showNorm=None, forceNorm=None, showMomentum=None, showPotential=None, - scalePsi=None, scaleMom=None, zoomMom=None, scalePot=None, + scalePsi=None, scaleMom=None, zoomMom=None, scalePot=None, psiRepresentation=None, momRepresentation=None, drawClassical=None, drawClassicalTrace=None, drawExpected=None, drawExpectedTrace=None, forceRecreate = False, customPlot = None): @@ -260,9 +323,11 @@ def reset_plot(self, width=None, height=None, if showNorm != None: self.showNorm = showNorm; recreate = True if forceNorm != None: self.forceNorm = forceNorm if showMomentum != None: self.showMomentum = showMomentum; recreate = True - if showPotential != None: self.showPotential = showPotential; + if showPotential != None: self.showPotential = showPotential; recreate = True if scalePsi != None: self.scalePsi = scalePsi if scaleMom != None: self.scaleMom = scaleMom + if psiRepresentation != None: self.psiRepresentation = psiRepresentation; recreate = True + if momRepresentation != None: self.momRepresentation = momRepresentation; recreate = True if zoomMom != None: self.zoomMom = zoomMom if scalePot != None: self.scalePot = scalePot if drawClassical != None: self.drawClassical = drawClassical @@ -288,6 +353,9 @@ def reset_plot(self, width=None, height=None, if cplot: if extraSubplots == 1: self.axCustom = self.fig.add_subplot(3,3,6) if not self.customPlotFull else self.fig.add_subplot(1,3,3) + elif extraSubplots == 2 and self.showMomentum: + self.axCustom = self.fig.add_subplot(3, 3, 3) + cplot = 0 else: self.axCustom = self.fig.add_subplot(4, 3, 3) else: self.axCustom = None @@ -309,16 +377,17 @@ def reset_plot(self, width=None, height=None, self.axMomentum = None # First drawing - self.QSystem.modSquared() + ###self.QSystem.modSquared() ### Not generic title = "Espai de posicions" if self.scalePsi: title += " (color reescalat)" self.axPsi.set_title(title) self.axPsi.set_xlabel("x ({})".format(self.unit_dist)) self.axPsi.set_ylabel("y ({})".format(self.unit_dist)) - self.datPsi = self.axPsi.imshow(self.QSystem.psiMod.T, origin='lower', + self.datPsi = plotComplex(self.QSystem.psi,self.axPsi,self.psiRepresentation,(self.QSystem.x0, self.QSystem.xf, self.QSystem.y0, self.QSystem.yf)) + """self.axPsi.imshow(self.QSystem.psiMod.T, origin='lower', extent=(self.QSystem.x0, self.QSystem.xf, self.QSystem.y0, self.QSystem.yf), aspect='equal', cmap="viridis", - interpolation='none') # , animated=True) # Doesn't do anything for imshow? + interpolation='none') # , animated=True) # Doesn't do anything for imshow?""" if self.showPotential: mathPhysics.set2DMatrix(self.QSystem.X, self.QSystem.Y, self.QSystem.potential, self.potentialMat, @@ -329,7 +398,7 @@ def reset_plot(self, width=None, height=None, interpolation='none') # , animated=True) # self.fig.colorbar(datPot, ax=ax1, label="Potencial: Força constant cap a baix") - if cplot: + if self.customPlot is not None: self.datCustom = {} # First plot self.customPlot[0](self.axCustom, self.datCustom, self.QSystem, self.units) @@ -352,10 +421,10 @@ def reset_plot(self, width=None, height=None, self.axEnergy.set_ylim(minValE-visibilityDelta, maxValE+visibilityDelta) self.axEnergy.set_ylabel("({})".format(self.unit_energy)) if not self.showNorm: self.axEnergy.set_xlabel("t ({})".format(self.unit_time)) - self.lineK, = self.axEnergy.plot(self.TList, self.KList, + self.lineK, = self.axEnergy.plot(self.TList, self.KList, 'red', label="K") # , animated=True) #Doesn't seem to speed up anything. Just complicates things - self.lineV, = self.axEnergy.plot(self.TList, self.VList, label="V") # , animated=True) - self.lineE, = self.axEnergy.plot(self.TList, self.EList, 'grey', label="E") # , animated=True) + self.lineV, = self.axEnergy.plot(self.TList, self.VList, 'blue', label="V") # , animated=True) + self.lineE, = self.axEnergy.plot(self.TList, self.EList, 'purple', label="E") # , animated=True) #Colors before were gray here and nothing for others self.axEnergy.grid() self.axEnergy.legend() @@ -532,12 +601,15 @@ def update(self, frame, onlyDraw=False): if self.debugTime: print("Time step (x{:d}): {:12.8f}".format(max(1, int(self.dtAnim / self.dtSim)) if self.stepsPerFrame == 0 else self.stepsPerFrame , time.time() - t0), " (s)", end=' <> ') - self.QSystem.modSquared() + #self.QSystem.modSquared() if self.psiCutoff != None: self.QSystem.psiMod[self.QSystem.psiMod < self.psiCutoff] = None - if self.scalePsi: self.datPsi.set_clim(vmax=np.max(self.QSystem.psiMod.T), vmin=0.) + # this is not optimal. It would be nice to have a numpy function that finds minimum and maximum on same loop # https://stackoverflow.com/questions/12200580/numpy-function-for-simultaneous-max-and-min + if self.scalePsi and self.psiRepresentation!='phase': self.rescalePsi()#self.datPsi.set_clim(vmax=np.max(self.datPsi.get_array()), + # vmin=0. if self.psiRepresentation == 'mod2' else np.min(self.datPsi.get_array())) + #self.QSystem.psiMod.T - self.datPsi.set_data(self.QSystem.psiMod.T) + self.datPsi.set_data(plotComplexData(self.QSystem.psi,self.psiRepresentation))#self.QSystem.psiMod.T) changes.append(self.datPsi) if self.showPotential: @@ -582,7 +654,7 @@ def update(self, frame, onlyDraw=False): if self.customPlot is not None: self.redraw = self.customPlot[1](self.axCustom, self.datCustom, self.QSystem, self.units) for dat in self.datCustom.values(): - changes.append(dat) + if issubclass(type(dat), matplotlib.artist.Artist): changes.append(dat) self.TList.append(self.QSystem.t) if self.showEnergy or self.forceEnergy: #Always store it? @@ -704,52 +776,53 @@ def manualUpdate(self, onlyDraw=False): def updatePotentialDraw(self, *args): # Updates potential but other stuff too. For example for visualizing changes in parameters #t0 = time.time() - self.updating = True - mathPhysics.set2DMatrix(self.QSystem.X, self.QSystem.Y, - self.QSystem.potential, self.potentialMat, - t=self.QSystem.t, extra_param=self.QSystem.extra_param) - if self.potentialCutoff != None: - self.potentialMat[self.potentialMat < self.potentialCutoff] = None - # if onlyDraw: self.datMom.set_clim(vmax=np.max(self.QSystem.psiMod.T), vmin=0.) # Being general we are slower! - if self.scalePot: self.datPot.set_clim(vmax=np.max(self.potentialMat), vmin=np.min(self.potentialMat)) # CAN BE OPTIMIZED - self.datPot.set_data(self.potentialMat.T) - - self.fig.draw_artist(self.datPsi) - self.fig.draw_artist(self.datPot) - - artists = set() # Using a set makes sure we don't repeat draw stuff - for patch in self.axPsi.patches: - #self.fig.draw_artist(patch) - artists.add(patch) - for line in self.axPsi.lines: - #self.fig.draw_artist(line) - artists.add(line) - for coll in self.axPsi.collections: - #self.fig.draw_artist(coll) - artists.add(coll) - - if self.extraUpdatesUpdate and self.extraUpdates is not None: - for action in self.extraUpdates: - updated = action(instance=self) - # except: updated = action() - if updated != None: - if type(updated) is tuple: - for el in updated: - artists.add(el) - else: - artists.add(updated) - if self.customPlotUpdate and self.customPlot is not None: - self.customPlot[1](self.axCustom, self.datCustom, self.QSystem, self.units) - for dat in self.datCustom.values(): - artists.add(dat) - - for art in artists: - self.fig.draw_artist(art) + if self.showPotential: + self.updating = True + mathPhysics.set2DMatrix(self.QSystem.X, self.QSystem.Y, + self.QSystem.potential, self.potentialMat, + t=self.QSystem.t, extra_param=self.QSystem.extra_param) + if self.potentialCutoff != None: + self.potentialMat[self.potentialMat < self.potentialCutoff] = None + # if onlyDraw: self.datMom.set_clim(vmax=np.max(self.QSystem.psiMod.T), vmin=0.) # Being general we are slower! + if self.scalePot: self.datPot.set_clim(vmax=np.max(self.potentialMat), vmin=np.min(self.potentialMat)) # CAN BE OPTIMIZED + self.datPot.set_data(self.potentialMat.T) + + self.fig.draw_artist(self.datPsi) + self.fig.draw_artist(self.datPot) - - self.fig.canvas.drawOptimized() - #print("Took {:12.8f} s, in FPS: {:12.8f}".format(time.time()-t0, 1./(time.time()-t0))) - self.updating = False + artists = set() # Using a set makes sure we don't repeat draw stuff + for patch in self.axPsi.patches: + #self.fig.draw_artist(patch) + artists.add(patch) + for line in self.axPsi.lines: + #self.fig.draw_artist(line) + artists.add(line) + for coll in self.axPsi.collections: + #self.fig.draw_artist(coll) + artists.add(coll) + + if self.extraUpdatesUpdate and self.extraUpdates is not None: + for action in self.extraUpdates: + updated = action(instance=self) + # except: updated = action() + if updated != None: + if type(updated) is tuple: + for el in updated: + artists.add(el) + else: + artists.add(updated) + if self.customPlotUpdate and self.customPlot is not None: + self.customPlot[1](self.axCustom, self.datCustom, self.QSystem, self.units) + for dat in self.datCustom.values(): + if issubclass(type(dat), matplotlib.artist.Artist): artists.add(dat) + + for art in artists: + self.fig.draw_artist(art) + + + self.fig.canvas.drawOptimized() + #print("Took {:12.8f} s, in FPS: {:12.8f}".format(time.time()-t0, 1./(time.time()-t0))) + self.updating = False def resetSystem(self, QSystem): self.QSystem = QSystem @@ -758,9 +831,16 @@ def resetSystem(self, QSystem): self.reset_plot(forceRecreate=True) def rescalePsi(self): - self.datPsi.set_clim(vmax=np.max(self.datPsi.get_array()), vmin=0.) + vmax = np.max(self.datPsi.get_array()) + if self.psiRepresentation == 'mod2': + self.datPsi.set_clim(vmax=np.max(vmax), vmin=0.) + else: + vmin = np.min(self.datPsi.get_array()) + maxVal = max(vmax, abs(vmin)) + self.datPsi.set_clim(vmax=maxVal, vmin=-maxVal) def saveAnimation(self, outputName, type="gif"): + #https://stackoverflow.com/questions/19646520/memory-usage-for-matplotlib-animation # Depends on current directory from pathlib import Path #cwd = Path.cwd() # Depends on how main.py is run. If run from terminal diff --git a/QuantumMovement/crankNicolson/crankNicolson2D.py b/QuantumMovement/crankNicolson/crankNicolson2D.py index a5cdcffe..14377241 100644 --- a/QuantumMovement/crankNicolson/crankNicolson2D.py +++ b/QuantumMovement/crankNicolson/crankNicolson2D.py @@ -23,6 +23,12 @@ Fortran to Python: https://numpy.org/doc/stable/f2py/ """ +"""# TODO +# General charged particle hamiltonian +# Also add spin to this, for example one color spin up other color spin down + +""" + # IMPORTS import types import numpy as np @@ -33,6 +39,7 @@ from scipy.linalg import solve_banded#, solveh_banded from scipy.fft import dstn from scipy.fft import idstn +from scipy.special import eval_hermite @@ -67,6 +74,15 @@ # ELECTRÓ: Si fem servir unitats estàndard: # - distancia (Å), temps (fs), energia (eV) +# +# També podem fer servir, si les coneixem, unitats atòmics (a.u.) +# hbar = 1, me = 1, e = 1, 4π€0 = 1 +# +# 1 d'energia : 1 Hartree ( Eh = 27.211 eV ) +# 1 de distància: 1 radi bohr ( 0.5292 Å ) +# 1 de temps : ( 0.0241888 fs ) +# +# @@ -74,7 +90,7 @@ # It's hard to test when is a method better than another # Sometimes abs2 is faster: Mod[:,:] = abs2(psi) -# Sometimes, when changing psi before, or with very big matrices, abs2Matrix is faster (jit functions have overhead). +# Sometimes, when changing psi before, or with very big matrices, abs2Matrix is faster (jit functions have overhead?). @numba.vectorize([numba.float64(numba.complex128), numba.float32(numba.complex64)]) def abs2(x): return x.real**2 + x.imag**2 @@ -407,7 +423,7 @@ def applyOperator2DOp(X, Y, psi, result, operator, doConjugate = False): if len(operator) != psi.ndim+1: raise ValueError("Operator doesn't match dimensions of matrix") if not doConjugate: for i in range(1, len(psi)-1): - for j in range(1, len(psi[0]-1)): + for j in range(1, len(psi[0])-1): result[i,j] = operator[0,0]*psi[i,j] + operator[1,0] * psi[i-1,j] + operator[1,1] * psi[i+1,j] + \ operator[2,0] * psi[i,j-1] + operator[2,1] * psi[i,j+1] result[0, :] = psi[0, :]*operator[0,0] @@ -416,7 +432,7 @@ def applyOperator2DOp(X, Y, psi, result, operator, doConjugate = False): result[:, len(psi[0])-1] = psi[:, len(psi[0])-1] * operator[0,0] else: for i in range(1, len(psi)-1): - for j in range(1, len(psi[0]-1)): + for j in range(1, len(psi[0])-1): result[i,j] = operator[0,0]*psi[i,j] + operator[1,0] * psi[i-1,j] + operator[1,1] * psi[i+1,j] + \ operator[2,0] * psi[i,j-1] + operator[2,1] * psi[i,j+1] result[i,j] = result[i,j] * np.conj(psi[i,j]) @@ -470,7 +486,7 @@ def applyOperator2DFuncNoJit(X, Y, psi, result, operator, t=0., extra_param=np.a if len(X) != len(psi) or len(Y) != len(psi[0]): raise ValueError("X and Y coordinates don't match shape of Psi") dx = (X[-1] - X[0]) / (len(psi) - 1) dy = (Y[-1] - Y[0]) / (len(psi[0]) - 1) - op = np.array([[1., 0.], [0., 0.], [0., 0.]]) + op = np.array([[1., 0.], [0., 0.], [0., 0.]], dtype=np.complex128) operator(op, (X[0], Y[0]), (dx, dy), t, extra_param=extra_param, onlyUpdate=False) if not doConjugate: for i in range(1, len(psi) - 1): @@ -906,6 +922,26 @@ def crankNicolson2DSchrodingerStepImaginary(X, Y, t, dt, potential, psi, psiTemp tridiagReal(gamma, tridY, tempBY, psi[itx,1:-1]) +@jit(nopython=True) +def eulerSchrodingerStep(X, Y, t, dt, potential, psi, psiTemp, extra_param=np.array([])): + """ Solves System: i h_red d/dt psi = [-h_red^2 /2 (d^2/dx^2 + d^2/dy^2) + V ] psi + Step 0: Psi(t) -> Step 1: psiTemp = Psi(t+dt/2) -> Step 2: psi = Psi(t+dt)""" + """ Same as VaryingPotential but Potential is taken only at t, as it's assumed constnat. Actually + to be used for imaginary time displacements, to find eigenstates. dt is imaginary component (real)""" + global hred, M + invdx2 = 1/( (X[-1] - X[0]) / (len(X) - 1) )**2 + invdy2 = 1/( (Y[-1] - Y[0]) / (len(Y) - 1) )**2 + for i in range(1, len(psi) - 1): + for j in range(1, len(psi[0]) - 1): + psi[i,j] += -dt * 1j/hred * ( -hred**2 / (2*M) * ( + (psi[i+1,j]-2*psi[i,j]+psi[i-1,j])*invdx2 + + (psi[i,j+1]-2*psi[i,j]+psi[i,j-1])*invdy2)+ + potential(X[i],Y[j],t, extra_param=extra_param)*psi[i,j]) + """psi[:,0] = 0. + psi[:,-1]= 0. + psi[0,:] = 0. + psi[-1,:]= 0.""" + @jit(nopython=True) def crankNicolson2DSchrodingerStepFastest(X, Y, t, dt, potential, psi, psiTemp, psiMod, extra_param=np.array([])): """ Solves System: i h_red d/dt psi = [-h_red^2 /2 (d^2/dx^2 + d^2/dy^2) + V ] psi @@ -1110,7 +1146,7 @@ def potential0(x, y, t=0., extra_param=np.array([])): def func1(x, y, t=0., extra_param=np.array([])): - return 1. + return 1#1/2. + np.sign(x) + np.sign(y) @jit def potentialBarrier(x, y, t, extra_param): @@ -1197,7 +1233,7 @@ def gaussian00(x,y,t=0., extra_param=np.array([])): def gaussianPacket(x, x0, sigma, p0, extra_param=np.array([])): global hred - return 1./(2*np.pi*sigma**2)**(0.25) * np.exp(-1./4. * ((x-x0)/sigma)**2) * np.exp(1j/hred * p0*(x))#-x0/2)) ??? Apunts mec quantica + return 1./(2*np.pi*sigma**2)**(0.25) * np.exp(-1./4. * ((x-x0)/sigma)**2) * np.exp(1j/hred * p0*(x))#-x0/2)) # Apunts mec quantica # GENERATOR FOR GAUSSIAN INITIAL STATES @@ -1217,16 +1253,17 @@ def inicial2D(x, y, t=0., extra_param=np.array([])): return inicial(x, p_0)*inicial(y, 0.5*p_0) -def eigenvectorsHarmonic1D(x, x0, n, k): +def eigenvectorsHarmonic1D(x, n, k): #not normalized x = x * np.sqrt(np.sqrt(M*k)/hred) - return np.power(M*k, 1/8.)/np.sqrt(2**n * np.math.factorial(n) * np.sqrt(hred*np.pi)) * np.exp(-(x-x0)**2/2.)*\ - np.polynomial.hermite.hermval(x, [0]*(n-1) + [1]) + return np.power(M*k, 1/8.)/np.sqrt(2**n * np.math.factorial(n) * np.sqrt(hred*np.pi)) * np.exp(-(x)**2/2.)*\ + eval_hermite(n, x * np.sqrt(np.sqrt(M*k)/hred)) + #np.polynomial.hermite.hermval(x, [0]*(n-1) + [1]) # GENERATOR FOR HARMONIC OSCILLATOR EIGENVECTORS def eigenvectorHarmonic2DGenerator(x0, nx, y0, ny, k): def result(x, y, t=0., extra_param=np.array([])): - return eigenvectorsHarmonic1D(x, x0, nx, k)*eigenvectorsHarmonic1D(y, y0, ny, k) + return eigenvectorsHarmonic1D(x-x0, nx, k)*eigenvectorsHarmonic1D(y-y0, ny, k) return result @@ -1346,7 +1383,7 @@ def evolveStep(self, dt): if self.step == 'eigen': crankNicolson2DSchrodingerStepClosedBoxEigen(self.X, self.Y, self.t, dt, self.potential, self.psi, self.psiEigen, self.psiMod, extra_param=self.extra_param) - if self.step == 'exact': + elif self.step == 'exact': crankNicolson2DSchrodingerStepExact(self.X, self.Y, self.t, dt, self.potential, self.psi, self.psiCopy, self.psiMod, extra_param=self.extra_param) else: @@ -1363,10 +1400,11 @@ def evolveImagStep(self, dt): """Potential will be taken as constant. Useful for imaginary time displacements, approximate eigenstates. NEGATIVE!""" crankNicolson2DSchrodingerStepImaginary(self.X, self.Y, self.t, dt, self.potential, self.psi, self.psiCopy, extra_param=self.extra_param) - + #eulerSchrodingerStep(self.X, self.Y, self.t, 1j*dt, self.potential, self.psi, self.psiCopy, extra_param=self.extra_param) # euler doesn't work # To avoid crashing in obtaining very very low values or very very high values # we normalize here - self.psi[:,:] = self.psi[:,:]/np.sqrt(self.norm()) + self.renorm() + #self.psi[:,:] = self.psi[:,:]/np.sqrt(self.norm()) def momentumSpace(self): fourierTransform2D(self.X, self.Y, self.psi, self.psiMom) @@ -1416,21 +1454,31 @@ def potentialEnergy(self): ####return np.trapz(np.trapz(self.psiMod))*self.dx*self.dy #return simpson_complex_list2D(self.dx, self.dy, self.psiMod) - def approximateEigenstate(self, maxiter = 100, callback=None, tol = 1e-4, resetInit=True): + def totalEnergy(self): + return self.kineticEnergy() + self.potentialEnergy() + + def approximateEigenstate(self, maxiter = 100, callback=None, tol = 1e-4, resetInit=True, initState=None): """ Iterates some time until, approximately, the newest lowest energy eigenstate remains Returns true if it was found u """ + """ + WE CAN CHANGE IT! WE DON'T NEED OT USE CRANK NICOLSON. MAYBE EULER METHOD FASTER AND CN NOT NECESSARY + Actually... Euler doesn't work + """ #print(len(self.eigenstates)) - - if resetInit: self.setState(func1) + if resetInit: + if initState is not None: + self.setState(initState) + else: + self.setState(func1) # PROBLEM WITH THIS. SYMMETRIC INITIAL STATE. HARD TO FIND ANTISYMMETRIC STATES for it in range(maxiter): for E, eigenstate in self.eigenstates: self.substractComponent(eigenstate) for __ in range(10): self.evolveImagStep(-2**(-6)) - self.renorm() + #self.renorm() # redundant, we put renormalization in evolveImagStep isEigen, E = self.isEigenstate(tol=tol) if isEigen: @@ -1737,6 +1785,7 @@ def crankNicolson2DHalfStepSchrodinger(X, Y, t, dt, potential, operator, psi, ps def aharonovBohmOperator(op, x, y, dx, dy, t=0, extra_param=np.array([1,1,1,1,0.2]), dir=-1, onlyUpdate=True): # -i/h [H-V] = -ih/2m [∂^2_x + ∂^2_y + i alpha / r^2 * (-y∂x + x∂y) - a^2/4r^2] alpha = extra_param[4] + x = x-1 #if onlyUpdate: return # The operator changes at every point invR2 = 1/(x*x+y*y+1e-10) diff --git a/QuantumMovement/main.py b/QuantumMovement/main.py index 3052f156..c96afc4c 100644 --- a/QuantumMovement/main.py +++ b/QuantumMovement/main.py @@ -1,25 +1,34 @@ +""" +TODO: +Aharonov Bohm / Double slit: Allow closing some of the slits +Sliders: Allow manual input of specific value. Ex 2.00, instead of 2.01, 1.99 +""" +from matplotlib import use -#import matplotlib -#matplotlib.use('agg') -#matplotlib.use('module://kivy.garden.matplotlib.backend_kivyagg') +use('Agg') # STOP MEMORY LEAKS #https://stackoverflow.com/questions/72271763/matplotlib-memory-leak-when-saving-figure-in-a-loop +# matplotlib.use('module://kivy.garden.matplotlib.backend_kivyagg') # Should this be used? Can't save to file then? Better not, implicitly draws figure using agg from functools import partial + # Combine functions. For example useful for callbacks def combF(f1, f2): def f(*args, **kwargs): f1(*args, **kwargs) f2(*args, **kwargs) + import random import json from pathlib import Path + mainDir = Path(__file__).parent from kivy.config import Config -Config.set('input', 'mouse', 'mouse,multitouch_on_demand') # Red circles appear when right clicking without this??? + +Config.set('input', 'mouse', 'mouse,multitouch_on_demand') # Red circles appear when right clicking without this??? import kivy from kivy.app import App from kivy.uix.label import Label @@ -43,13 +52,12 @@ def f(*args, **kwargs): from kivy.metrics import dp, sp - from kivy.properties import ObjectProperty from kivy.properties import NumericProperty from kivy.garden.matplotlib.backend_kivyagg import FigureCanvasKivyAgg -#from kivy.garden.matplotlib.backend_kivyagg import FigureCanvas +# from kivy.garden.matplotlib.backend_kivyagg import FigureCanvas from crankNicolson.animate import FigureCanvasKivyAggModified -#from kivy.lang import Builder +# from kivy.lang import Builder from kivy.uix.screenmanager import ScreenManager, Screen from kivy.clock import Clock from kivy.factory import Factory @@ -60,22 +68,25 @@ def f(*args, **kwargs): import crankNicolson.crankNicolson2D as mathPhysics import crankNicolson.animate as animate import matplotlib.pyplot as plt + plt.rcParams.update({'font.size': 11}) #### # https://stackoverflow.com/questions/70629758/kivy-how-to-pass-arguments-to-a-widget-class ##### -#import warnings # For debugging. Warnings stop the program, and thus also show where they occur -#warnings.filterwarnings('error') -#import matplotlib.cbook.deprecation -#warnings.filterwarnings('ignore', category=matplotlib.cbook.deprecation.MatplotlibDeprecationWarning) +# import warnings # For debugging. Warnings stop the program, and thus also show where they occur +# warnings.filterwarnings('error') +# import matplotlib.cbook.deprecation +# warnings.filterwarnings('ignore', category=matplotlib.cbook.deprecation.MatplotlibDeprecationWarning) +# import gc + +unit_dist = r"$a_0$" # '2 Å'#r"a_0"# +unit_time = r"$\frac{\hbar}{E_h}$" # '1/3 fs'#r""# +unit_energy = r"$E_h$" # '2 eV' #r"E_h"# +unit_mom = r"$\frac{\hbar}{a_0}$" # r'$\frac{1}{2}\hbar Å^{-1}$'#'1/3 eV·fs/Å' #'2 eV · 1/3 fs / 2 Å' -unit_dist = '2 Å' -unit_time = '1/3 fs' -unit_energy = '2 eV' -unit_mom = r'$\frac{1}{2}\hbar Å^{-1}$'#'1/3 eV·fs/Å' #'2 eV · 1/3 fs / 2 Å' # hred = 2 eV · 1/3 fs = 2/3 eV · fs # ℏ ℏ ℏ ℏ @@ -83,73 +94,95 @@ def f(*args, **kwargs): # We could be more general and just define a general button with image, # but we use them for just particular cases anyway +class LightButton(Button): + pass + + +class LightButtonImage(LightButton): + def __init__(self, image_src, size_hint_image=1., **kwargs): + self.image_src = image_src + self.size_hint_image = size_hint_image + super(LightButtonImage, self).__init__(**kwargs) + + pass + + class PlayButton(ToggleButton): pass + class SaveButton(Button): pass + class ReturnButton(Button): pass + class SettingsButton(Button): pass + class InfoButton(Button): pass + class ColoredGridLayout(GridLayout): pass + class ColoredBoxLayout(BoxLayout): pass + class BoolCheckBox(CheckBox): pass + class FastGraphicsCheckbox(CheckBox): def on_active(self, *args): animate.optimizeGraphics = self.active - - class GlobalVariable(GridLayout): def __init__(self, source, names, sliders, num=0, **kwargs): - super(GlobalVariable, self).__init__(cols=4,**kwargs) + super(GlobalVariable, self).__init__(cols=4, **kwargs) self.num = num self.names = names - self.sliders = sliders # List of: show variable i as a slider? - self.source = source # List of variables, numpy array + self.sliders = sliders # List of: show variable i as a slider? + self.source = source # List of variables, numpy array self.label = Label(text="Var.\n{}".format(num), size_hint_x=0.2) - self.nameText = DataInput(attribute="names", index=self.num, holder=self, condition="unique", multiline=False, size_hint_x=0.3) + self.nameText = DataInput(attribute="names", index=self.num, holder=self, condition="unique", multiline=False, + size_hint_x=0.3) self.valText = DataInput(attribute="source", index=self.num, holder=self, multiline=False, size_hint_x=0.3) self.sliderQuery = GridLayout(rows=2, size_hint_x=0.2) self.sliderQuery.add_widget(Label(text="Slider?")) sliderCheck = CheckBox(active=sliders[num]) + def updateSlider(checkbox, active): sliders[num] = active + sliderCheck.bind(active=updateSlider) self.sliderQuery.add_widget(sliderCheck) - - - #self.layout = GridLayout(cols=3, size=(100, 100)) + # self.layout = GridLayout(cols=3, size=(100, 100)) self.add_widget(self.label) self.add_widget(self.nameText) self.add_widget(self.valText) self.add_widget(self.sliderQuery) - #self.add_widget(self.layout) + # self.add_widget(self.layout) """def on_value(self, instance, value): self.source[self.num] = self.value""" + class GlobalVariablesPopup(Popup): nVar = 16 + def __init__(self, window, **kwargs): self.window = window # Window holds information such as QuantumSystem and Animation super(GlobalVariablesPopup, self).__init__(**kwargs) @@ -158,17 +191,20 @@ def __init__(self, window, **kwargs): Clock.schedule_once(self._finish_init) def _finish_init(self, *args): - #print(self.window.paramNames) + # print(self.window.paramNames) for i in range(self.nVar): - self.layout.add_widget(GlobalVariable(self.window.extra_param, self.window.paramNames, self.window.paramSliders, num=i)) + self.layout.add_widget( + GlobalVariable(self.window.extra_param, self.window.paramNames, self.window.paramSliders, num=i)) self.add_widget(self.layout) def on_dismiss(self): super(GlobalVariablesPopup, self).on_dismiss() + # We need to wait to make sure all DataInputs finish on unfocus def updateStuff(*param): self.window.setVarDict() self.window.setSliders() + Clock.schedule_once(updateStuff) @@ -183,24 +219,22 @@ def _finish_init(self, *args): # Lambda in loops! Careful with using iterating variable # Same with function, using global variable instead of value during creation # https://stackoverflow.com/questions/19837486/lambda-in-a-loop - #grid = self.ids.states - #self.parts = [] + # grid = self.ids.states + # self.parts = [] for state in self.window.savedStates: self.add_state(state) - - # Like this we can test everything gets correctly garbage collected (eventually) def on_dismiss(self): super(SavedStatesPopup, self).on_dismiss() - #gc.collect() + # gc.collect() def add_state(self, state): grid = self.ids.states lbl = Label(text=state["name"]) # , btnprev = Button(text="Previsualitza", - on_release=lambda x, state=state: PlotPopup(state).open()) + on_release=lambda x, state=state: Factory.PlotPopup(state).open()) btnchan = Button(text="Canvia\na aquest", on_release=lambda x, state=state: self.window.setState(state)) btnsub = Button(text="Substreu\na l'actual", @@ -221,7 +255,7 @@ def removeBind(*args, state=state, grid.remove_widget(btnchan) grid.remove_widget(btnsub) grid.remove_widget(btndel) - self.window.ids.stateName.text = "estat{}".format(len(self.window.savedStates)) + self.window.ids.stateName.text = "est{}".format(len(self.window.savedStates)) btndel.bind(on_release=removeBind) @@ -231,15 +265,58 @@ def removeBind(*args, state=state, grid.add_widget(btnsub) grid.add_widget(btndel) + class SavedEigenstatesPopup(Popup): def __init__(self, window, **kwargs): self.window = window - self.tol = 1e-4 + self.tol = 1e-6 self.maxiter = 20 + self.isClose = False super(SavedEigenstatesPopup, self).__init__(**kwargs) Clock.schedule_once(self._finish_init) + def add_state(self, state, grid, E, count=0, index=0): + lbl = Label(text="E={:.4e} |{}".format(E, count)) # , + + stateCopy = state.copy() + btnprev = Button(text="Previsualitza", + on_release=lambda x, state=stateCopy: Factory.PlotPopup(state).open()) + btnchan = Button(text="Canvia\na aquest", + on_release=lambda x, state=stateCopy: self.window.setState(state)) + btnsub = Button(text="Substreu\na l'actual", + on_release=lambda x, state=stateCopy: self.window.substractComponent(state)) + + # Bad, we are multiplying memory use here... We should not create more copies of states for each button + + btndel = Button(text="Elimina", background_color=(0.6, 0, 0, 0.8)) + + def removeBind(*args, state=state, + lbl=lbl, btnprev=btnprev, btnchan=btnchan, btnsub=btnsub, btndel=btndel, E=E, count=count): + # print(state["name"]) + rep = 0 + for j in range(len(self.window.QSystem.eigenstates)): + if self.window.QSystem.eigenstates[j][0] == E: + if rep == count: + del self.window.QSystem.eigenstates[j] + break + else: + rep += 1 + + grid.remove_widget(lbl) + grid.remove_widget(btnprev) + grid.remove_widget(btnchan) + grid.remove_widget(btnsub) + grid.remove_widget(btndel) + + btndel.bind(on_release=removeBind) + + grid.add_widget(lbl, index) + grid.add_widget(btnprev, index) + grid.add_widget(btnchan, index) + grid.add_widget(btnsub, index) + grid.add_widget(btndel, index) + def _finish_init(self, *args): """def updateProgress(val): @@ -250,7 +327,7 @@ def _finish_init(self, *args): # Same with function, using global variable instead of value during creation # https://stackoverflow.com/questions/19837486/lambda-in-a-loop grid = self.ids.states - #self.parts = [] + # self.parts = [] count = 0 prev = None @@ -258,13 +335,13 @@ def _finish_init(self, *args): self.window.tempState["psi"] = eigen state = self.window.tempState - if prev == E: count+=1 - lbl = Label(text="E={:.4e} |{}".format(E, count))#, + if prev == E: count += 1 + lbl = Label(text="E={:.4e} |{}".format(E, count)) # , prev = E stateCopy = state.copy() btnprev = Button(text="Previsualitza", - on_release = lambda x, state=stateCopy: PlotPopup(state).open()) + on_release=lambda x, state=stateCopy: Factory.PlotPopup(state).open()) btnchan = Button(text="Canvia\na aquest", on_release=lambda x, state=stateCopy: self.window.setState(state)) btnsub = Button(text="Substreu\na l'actual", @@ -272,18 +349,19 @@ def _finish_init(self, *args): # Bad, we are multiplying memory use here... We should not create more copies of states for each button - btndel = Button(text="Elimina", background_color=(0.6,0,0,0.8)) + btndel = Button(text="Elimina", background_color=(0.6, 0, 0, 0.8)) def removeBind(*args, state=state, lbl=lbl, btnprev=btnprev, btnchan=btnchan, btnsub=btnsub, btndel=btndel, E=E, count=count): - #print(state["name"]) + # print(state["name"]) rep = 0 for j in range(len(self.window.QSystem.eigenstates)): - if self.window.QSystem.eigenstates[j][0]==E: + if self.window.QSystem.eigenstates[j][0] == E: if rep == count: del self.window.QSystem.eigenstates[j] break - else: rep+=1 + else: + rep += 1 grid.remove_widget(lbl) grid.remove_widget(btnprev) @@ -298,28 +376,43 @@ def removeBind(*args, state=state, grid.add_widget(btnchan) grid.add_widget(btnsub) grid.add_widget(btndel) + def eigenFind(self): self.count = 0 - self.window.QSystem.setState(mathPhysics.func1) + if not self.isClose: self.window.QSystem.setState(mathPhysics.func1) self.window.QSystem.renorm() + def eigenLoop(*args): if self.count < self.maxiter: if self.window.QSystem.approximateEigenstate(tol=self.tol, maxiter=5, resetInit=False): self.window.animation.manualUpdate(onlyDraw=True) - self.dismiss() + # self.dismiss() + qsys = self.window.QSystem + E = np.real(qsys.expectedValueOp(qsys.totalEnergyOp, + doConjugate=False)) # np.real(self.window.QSystem.totalEnergy()) + count = -1 + index = 0 + for Et, eigen in self.window.QSystem.eigenstates[::-1]: + if Et > E: index += 5 + if Et == E: count += 1 # one time it is always found + # print(index) + self.window.tempState["psi"] = np.copy(qsys.psi) + self.add_state(self.window.tempState, self.ids.states, E, count=count, index=index) + self.isClose = False # We already found state. Now we are far form other eigenstates (orthogonal) + self.ids.eigenButton.text = "Busca següent\nEigenstat" return else: Clock.schedule_once(eigenLoop, 0.) - self.count+=1 - self.ids.progress.value = self.count/self.maxiter * 100 + self.count += 1 + self.ids.progress.value = self.count / self.maxiter * 100 else: TextPopup("No s'ha pogut trobar!").open() self.window.animation.manualUpdate(onlyDraw=True) - self.dismiss() - Clock.schedule_once(eigenLoop) - - + self.isClose = True + self.ids.eigenButton.text = "Continua\nbuscant" + # self.dismiss() + Clock.schedule_once(eigenLoop) class DataInput(TextInput): @@ -330,7 +423,9 @@ class DataInput(TextInput): Which returns a dictionary of variables, mutable. """ text_width = NumericProperty() - def __init__(self, attribute = None, index=None, holder = None, condition = None, callback=None, centered=False, scientific=False, maxDecimal=6, **kwargs): + + def __init__(self, attribute=None, index=None, holder=None, condition=None, callback=None, centered=False, + scientific=False, maxDecimal=6, **kwargs): self.scientific = False self.maxDecimal = maxDecimal self.centered = centered @@ -348,11 +443,11 @@ def __init__(self, attribute = None, index=None, holder = None, condition = None def set_self(self, dt): # Attribute is set in kivy language - #print(id(self.attribute)) Check - self.attributeVal = getattr(self.holder,self.attribute) if self.index is None else \ - getattr(self.holder,self.attribute)[self.index] + # print(id(self.attribute)) Check + self.attributeVal = getattr(self.holder, self.attribute) if self.index is None else \ + getattr(self.holder, self.attribute)[self.index] self.type = type(self.attributeVal) - form = "{:." + str(self.maxDecimal) + ('f' if not self.scientific else 'e') +'}' + form = "{:." + str(self.maxDecimal) + ('f' if not self.scientific else 'e') + '}' self.text = str(self.attributeVal) if self.type is not float else form.format(self.attributeVal).rstrip('0') self.copy = self.text @@ -362,21 +457,19 @@ def _on_focus(self, instance, value, *largs): if value == False and self.text != self.copy: self.on_text_validate() - def update_padding(self, *args): ''' Update the padding so the text is centered ''' self.text_width = self._get_text_width(self.text, self.tab_width, self._label_cached) - def on_text_validate(self): # Enter try: self.attributeVal = self.type(self.text) # Try first to see if it can be converted sucessfully except: - #print("Couldn't convert properly") + # print("Couldn't convert properly") TextPopup("Invalid Input!").open() self.text = self.copy return @@ -392,7 +485,8 @@ def on_text_validate(self): self.callback(self.attributeVal) def conditionHolds(self, val): - if self.condition == None: return True + if self.condition == None: + return True elif self.condition == "fixed": TextPopup("Fixed Value", title="Warning").open() return False @@ -401,7 +495,7 @@ def conditionHolds(self, val): elif self.condition == "unique" and self.index != None: ### Not too general ### We assume here no name is allowed to be repeated as a name - if val != "" and val in getattr(self.holder,self.attribute): + if val != "" and val in getattr(self.holder, self.attribute): TextPopup("Can't repeat!", title="Warning").open() return False elif self.condition == "nonnegative": @@ -417,37 +511,36 @@ def conditionHolds(self, val): left = self.condition.split('-')[1] right = self.condition.split('-')[2] if not (float(left) <= val <= float(right)): - #if val < float(left) or float(right) < val: + # if val < float(left) or float(right) < val: TextPopup("Must be between {0} and {1}".format(left, right), title="Warning").open() return False elif self.condition.startswith("gt"): if not (float(self.condition[2:]) < val): - TextPopup("Must be greater than "+self.condition[2:], title="Warning").open() + TextPopup("Must be greater than " + self.condition[2:], title="Warning").open() return False elif self.condition.startswith("geq"): if not (float(self.condition[2:]) <= val): - TextPopup("Must be greater than or equal to "+self.condition[2:], title="Warning").open() + TextPopup("Must be greater than or equal to " + self.condition[2:], title="Warning").open() return False elif self.condition.startswith("lt"): if not (val < float(self.condition[2:])): - TextPopup("Must be less than "+self.condition[2:], title="Warning").open() + TextPopup("Must be less than " + self.condition[2:], title="Warning").open() return False elif self.condition.startswith("leq"): if not (val <= float(self.condition[2:])): - TextPopup("Must be less than or equal to "+self.condition[2:], title="Warning").open() + TextPopup("Must be less than or equal to " + self.condition[2:], title="Warning").open() return False return True - - class DataSlider(Slider): """ Similar and linked to a DataInput, but user fiendly form of a Slider. Only numeric variables! """ + # By default step = 0, which means pixel resolution for changes - def __init__(self, attribute = None, index=None, holder = None, callback=None, isPotential=True, **kwargs): + def __init__(self, attribute=None, index=None, holder=None, callback=None, isPotential=True, **kwargs): self.attribute = attribute self.index = index self.holder = holder @@ -463,7 +556,7 @@ def set_self(self, dt): # Attribute is set in kivy language # print(id( self.attribute)) Check self.attributeVal = getattr(self.holder, self.attribute) if self.index is None else \ - getattr(self.holder,self.attribute)[self.index] + getattr(self.holder, self.attribute)[self.index] self.text = str(self.attributeVal) self.copy = self.text @@ -478,9 +571,12 @@ def on_value(self, instance, val): self.callback(val) if not self.firstTime: - if self.isPotential and self.holder.animation.paused or self.holder.animation.isOver: self.holder.animation.updatePotentialDraw() - else: self.holder.animation.updating = True - else: self.firstTime = not self.firstTime + if self.isPotential and self.holder.animation.paused or self.holder.animation.isOver: + self.holder.animation.updatePotentialDraw() + else: + self.holder.animation.updating = True + else: + self.firstTime = not self.firstTime class CustomDataSlider(BoxLayout): @@ -488,7 +584,7 @@ class CustomDataSlider(BoxLayout): # Slider + custom slider ranges (max val and min val) ###### Maybe Clock once???? Need to check, maybe clock_once already inside slider/datainput is enough - def __init__(self, name = None, attribute = None, index=None, holder = None, orientation="horizontal", min = -1., max = 1., + def __init__(self, name=None, attribute=None, index=None, holder=None, orientation="horizontal", min=-1., max=1., isPotential=True, value=None, step=0., variableLimits=True, callback=None, **kwargs): super(CustomDataSlider, self).__init__(orientation=orientation, **kwargs) @@ -501,8 +597,8 @@ def __init__(self, name = None, attribute = None, index=None, holder = None, ori self.index = index self.value = value self.min = min - self.max=max - self.orientation=orientation + self.max = max + self.orientation = orientation self.isPotential = isPotential self.step = step self.callback = callback @@ -513,15 +609,17 @@ def _finish_init(self, dt): if isHor: self.size_hint_x = 1 self.size_hint_y = None - self.height = sp(80/2) + self.height = sp(80 / 2) sizeHintSlid = (0.7, 1) else: self.size_hint_x = None self.size_hint_y = 1 - self.width = sp(80/2) + self.width = sp(80 / 2) sizeHintSlid = (1, 1) - val = float(getattr(self.holder, self.attribute) if self.index is None else getattr(self.holder, self.attribute)[self.index]) + val = float( + getattr(self.holder, self.attribute) if self.index is None else getattr(self.holder, self.attribute)[ + self.index]) if self.value is not None: setVal = self.value @@ -539,11 +637,10 @@ def callback(newVal): if self.callback is not None: self.callback(newVal) - - self.slider = DataSlider(self.attribute, self.index, self.holder, value=setVal, min=self.min, max=self.max, size_hint=sizeHintSlid, - orientation=self.orientation, callback=callback, isPotential=self.isPotential, step=self.step) + orientation=self.orientation, callback=callback, isPotential=self.isPotential, + step=self.step) sizeHintText = (0.1, 1) if isHor else (1, None) self.minDat = DataInput(attribute="min", holder=self.slider, condition="lt{0}".format(self.slider.max), @@ -552,7 +649,6 @@ def callback(newVal): def updateMin(newMax): self.minDat.condition = "lt{0}".format(newMax) - self.maxDat = DataInput(attribute="max", holder=self.slider, condition="gt{0}".format(self.slider.min), size_hint=sizeHintText, centered=True, disabled=not self.variableLimits, height=sp(20)) @@ -563,7 +659,7 @@ def updateMax(newMin): self.minDat.callback = updateMax # Layout [min] [----- slider -----] [max] - #self.layout = BoxLayout(orientation=self.orientation) + # self.layout = BoxLayout(orientation=self.orientation) self.add_widget(self.label) if isHor: self.add_widget(self.minDat) @@ -573,11 +669,55 @@ def updateMax(newMin): self.add_widget(self.maxDat) self.add_widget(self.slider) self.add_widget(self.minDat) - #self.add_widget(self.layout) + # self.add_widget(self.layout) + + +################################################################################################################ +# --------------------------------------------------------------------------------------------------------------# +# FUNCIONS EN GENERAL # + +from crankNicolson.crankNicolson2D import gaussianPacket +from crankNicolson.crankNicolson2D import eigenvectorsHarmonic1D as eigenHarmonic + +""" +def gaussianPacket(x, x0, sigma, p0, extra_param=None): + global hred + return 1./(2*np.pi*sigma**2)**(0.25) * np.exp(-1./4. * ((x-x0)/sigma)**2) * np.exp(1j/hred * p0*(x)) +""" +@jit # (cache=True) +def heaviside(x, k=1.): + """ + Heaviside. Analytic approximation: 1/(1 + e^-2kr). Larger k, better approximation, less smooth + """ + return 1. / (1. + np.exp(-2. * k * x)) + + +# IMPORTANT +# IMPORTANT +# WARNING: eval() https://nedbatchelder.com/blog/201206/eval_really_is_dangerous.html +# WARNING +# This creates a GLOBAL function. This is why it's allowed to return +# it seems locals can't change, but globals can +# https://stackoverflow.com/questions/41100196/exec-not-working-inside-function-python3-x +from numpy import sin, cos, tan, arcsin, arccos, arctan, hypot, arctan2, degrees, radians, deg2rad, rad2deg, \ + sinh, cosh, tanh, arcsinh, arccosh, arctanh, \ + around, rint, fix, floor, ceil, trunc, \ + exp, expm1, exp2, log, log10, log2, log1p, logaddexp, logaddexp2, \ + power, sqrt, i0, sinc, \ + sign, \ + pi + +ln = log # for commodity +i = 1j + +# Heaviside is not supported in numba + sandbox = None + + class LocalOperator: """ A local operator is an operator that only acts on a point and its closest neighbours @@ -594,7 +734,8 @@ def __init__(self, operator): self.operator = operator else: def operatorFunc(op, X, dx, t=0, extra_param=np.array([]), dir=-1, onlyUpdate=True, operator=operator): - op[:,:] = operator[:,:] + op[:, :] = operator[:, :] + self.operator = operatorFunc # L * psi. Apply Operator @@ -602,35 +743,60 @@ def __mul__(self, state): state = state.state newState = state.copy() newState["psi"] = state["psi"].copy() - #print(newState["psi"]) + # print(newState["psi"]) sb = sandbox X = np.linspace(state["x0"], state["xf"], len(state["psi"])) Y = np.linspace(state["y0"], state["yf"], len(state["psi"][0])) - mathPhysics.applyOperator2DFuncNoJit(X, Y, state["psi"], newState["psi"], self.operator, t=sb.QSystem.t, extra_param=sb.QSystem.extra_param) + mathPhysics.applyOperator2DFuncNoJit(X, Y, state["psi"], newState["psi"], self.operator, t=sb.QSystem.t, + extra_param=sb.QSystem.extra_param) return OperableState(newState) # k L. Multiply by number def __rmul__(self, other): # numero - #if not callable(other): + # if not callable(other): + # sb = sandbox def operator(op, X, dx, t=0, extra_param=np.array([]), dir=-1, onlyUpdate=True): self.operator(op, X, dx, t, extra_param, dir, onlyUpdate) op *= other return LocalOperator(operator) + # Multiply by ndarray + def __array_ufunc__(self, ufunc, method, *inputs, **kwargs): + ndarray, loc_op = inputs + sb = sandbox + + def operator(op, X, dx, t=0, extra_param=np.array([]), dir=-1, onlyUpdate=True): + loc_op.operator(op, X, dx, t, extra_param, dir, onlyUpdate) + ix = round((X[0] - sb.QSystem.X[0]) / dx[0]) + iy = round((X[1] - sb.QSystem.Y[0]) / dx[1]) + op *= ndarray[ix, iy] + + return LocalOperator(operator) + # Add two operators. When we combine them we dont know if they overlap. WE'll havae to always update. # These local operators are only to be used to manipulate states. One operation at a time, not like kinetic energy def __add__(self, other): + if type(other) is np.ndarray: + other = other * LocalOperator(np.ndarray([[1., 0], [0, 0], [0, 0]])) # Other * Identity + def operator(op, X, dx, t=0, extra_param=np.array([]), dir=-1, onlyUpdate=False): self.operator(op, X, dx, t, extra_param, dir, onlyUpdate=False) opCopy = op.copy() other.operator(opCopy, X, dx, t, extra_param, dir, onlyUpdate=False) - op+=opCopy + op += opCopy + return LocalOperator(operator) + __radd__ = __add__ + + def __sub__(self, other): + return self.__add__(-1 * other) + + class OperableState(): # Dictionary holding state info def __init__(self, state): @@ -645,17 +811,36 @@ def __add__(self, other): def __mul__(self, other): newState = self.state.copy() newState["psi"] = self.state["psi"].copy() - newState["psi"]*=other + newState["psi"] *= other return OperableState(newState) + __rmul__ = __mul__ - #To multiply with ndarrays + def __truediv__(self, other): + return 1 / other * self + + def __matmul__(self, other): + """ + # SHOULD ONLY BE USED BETWEEN OPERABLE STATES! + # psi @ means using psi as a ket, + # ie. psi @ bra will evaluate + """ + dx = (self.state["xf"] - self.state["x0"]) / (len(self.state["psi"]) - 1) + dy = (self.state["yf"] - self.state["y0"]) / (len(self.state["psi"][0]) - 1) + return mathPhysics.innerProduct2D(self.state["psi"], other.state["psi"], dx, dy) + + def __sub__(self, other): + return self.__add__(-1 * other) + + # To multiply with ndarrays def __array_ufunc__(self, ufunc, method, *inputs, **kwargs): ndarray, state = inputs return state * ndarray + class StateExpression(TextInput): """Expected to hold string python expression which is to be evaluated.""" + rememberExpression = '' def __init__(self, varDict={}, holder=None, **kwargs): super(StateExpression, self).__init__(**kwargs) @@ -664,6 +849,8 @@ def __init__(self, varDict={}, holder=None, **kwargs): self.stateDict = {} + self.text = StateExpression.rememberExpression + Clock.schedule_once(self.set_self) # .kv can't take data from another class during init. Everything needs to be init first # That's why the delay @@ -671,42 +858,64 @@ def __init__(self, varDict={}, holder=None, **kwargs): def set_self(self, *args): self.states = self.holder.window.savedStates + def evaluate(self): + pass + def on_text_validate(self): # Enter + tempName = self.holder.ids.stateName.text + if tempName == "" or tempName[0] == ':' or tempName[0].isnumeric(): + Factory.TextPopup("Nom ha de començar amb\nuna lletra").open() + return + if not tempName.isidentifier(): + TextPopup("Nom invàlid!").open() + return + + for state in self.states: + if self.holder.ids.stateName.text == state["name"]: + TextPopup("Nom ja fet servir!").open() + return + + extra_param = self.holder.window.extra_param + self.stateDict.clear() stateVar = {} sbox = self.holder.window for state in self.states: + # Should we interpolate all states so we can operate on all of them? + # But cost of interpolating every single time!!! + # if state.get("x0") stateVar[state["name"]] = OperableState(state) - self.stateDict[state["name"]] = "stateVar['"+state["name"]+"']" - - # we assume we are working with current + self.stateDict[state["name"]] = "stateVar['" + state["name"] + "']" + # we assume we are working with current size Nx/Ny + # States that have been saved differently are interpolated K = LocalOperator(mathPhysics.kineticEnergy) - Px = LocalOperator(-1j * mathPhysics.hred * np.array([[0., 0.], [-1./(2.*sbox.QSystem.dx), +1./(2.*sbox.QSystem.dx)] , [0.,0.]]) ) - Dx = LocalOperator( np.array([[0., 0.], [-1./(2.*sbox.QSystem.dx), +1./(2.*sbox.QSystem.dx)] , [0.,0.]]) ) - Py = LocalOperator(-1j * mathPhysics.hred * np.array([[0., 0.], [0.,0.], [-1./(2.*sbox.QSystem.dy), +1./(2.*sbox.QSystem.dy)]])) - Dy = LocalOperator( np.array([[0., 0.], [0.,0.], [-1./(2.*sbox.QSystem.dy), +1./(2.*sbox.QSystem.dy)]])) + Px = LocalOperator(-1j * mathPhysics.hred * np.array( + [[0., 0.], [-1. / (2. * sbox.QSystem.dx), +1. / (2. * sbox.QSystem.dx)], [0., 0.]])) + Dx = LocalOperator(np.array([[0., 0.], [-1. / (2. * sbox.QSystem.dx), +1. / (2. * sbox.QSystem.dx)], [0., 0.]])) + Py = LocalOperator(-1j * mathPhysics.hred * np.array( + [[0., 0.], [0., 0.], [-1. / (2. * sbox.QSystem.dy), +1. / (2. * sbox.QSystem.dy)]])) + Dy = LocalOperator(np.array([[0., 0.], [0., 0.], [-1. / (2. * sbox.QSystem.dy), +1. / (2. * sbox.QSystem.dy)]])) + I = LocalOperator(np.array([[1., 0.], [0., 0.], [0., 0.]])) + H = LocalOperator(sbox.QSystem.totalEnergyOp) x = sbox.QSystem.Xmesh y = sbox.QSystem.Ymesh - repeated = False - for state in self.states: - if self.holder.ids.stateName.text == state["name"]: - repeated = True - break + L = (x * Py - y * Px) - if repeated: - TextPopup("Nom ja fet servir!").open() - return + V = np.copy(sbox.QSystem.psiMod) + mathPhysics.set2DMatrix(sbox.QSystem.X, sbox.QSystem.Y, + sbox.QSystem.potential, V, + t=sbox.QSystem.t, extra_param=sbox.QSystem.extra_param) expression = self.text if expression == "": return try: # Things like {px} are substituted by their corresponding actual extra_param - expressionFormated = expression.format(**self.stateDict) + expressionFormated = expression.format(**self.stateDict, **self.varDict) except: TextPopup("Compte amb estats o\nvariables globals").open() return @@ -719,25 +928,41 @@ def on_text_validate(self): # Molt malament try: - tempState = eval(expressionFormated).state.copy() - - tempState["name"] = self.holder.ids.stateName.text - self.states.append(tempState) - self.holder.add_state(tempState) + i = 1j + evaluated = eval(expressionFormated) + # print(evaluated) #debug + # print(type(evaluated)) + if type(evaluated) is not OperableState: # (complex, float, int): + TextPopup("Resultat numèric:\n{}".format(evaluated)).open() + else: + tempState = evaluated.state.copy() + + tempState["name"] = self.holder.ids.stateName.text + tempState["x0"] = sbox.QSystem.x0 + tempState["xf"] = sbox.QSystem.xf + tempState["y0"] = sbox.QSystem.y0 + tempState["yf"] = sbox.QSystem.yf + self.states.append(tempState) + self.holder.add_state(tempState) except: TextPopup("La declaració ha de\nresultar en un estat").open() - return 0 # Returns True when everything is OK + return 0 # Returns True when everything is OK def conditionHolds(self, val): if self.condition == None: return True return True + def on_text(self, object, value): # object will be self, it's quite redundant + # super(StateExpression, self).on_text(value, *args, **kwargs) #there is no default on text + StateExpression.rememberExpression = value + class FunctionInput(TextInput): """Expected to hold string python expression which can be converted into a function.""" - def __init__(self, functionName=None, definitionName=None, varDict={}, holder=None, condition=None, jit = False, **kwargs): + def __init__(self, functionName=None, definitionName=None, varDict={}, holder=None, condition=None, jit=False, + **kwargs): super(FunctionInput, self).__init__(**kwargs) self.functionName = functionName self.definitionName = definitionName @@ -769,86 +994,45 @@ def on_text_validate(self): setattr(self.holder, self.functionName, jit(self.func) if self.jit else self.func) setattr(self.holder, self.definitionName, self.definition) - return 0 # Returns True when everything is OK + return 0 # Returns True when everything is OK except MaliciousInput: exit(print("ALERTA: AIXÒ NO ES POT FER")) except InvalidFormat: TextPopup("Compte amb les variables globals").open() except: - TextPopup("Expressió Invàlida!\nRecorda multiplicar amb *\nI compte amb divisions per 0").open() #\nPer 1/r posa 1/(r+numPetit) - + TextPopup( + "Expressió Invàlida!\nRecorda multiplicar amb *\nI compte amb divisions per 0").open() # \nPer 1/r posa 1/(r+numPetit) def conditionHolds(self, val): if self.condition == None: return True return True -################################################################################################################ -#--------------------------------------------------------------------------------------------------------------# -# FUNCIONS EN GENERAL # - -from crankNicolson.crankNicolson2D import gaussianPacket -from crankNicolson.crankNicolson2D import eigenvectorsHarmonic1D as eigenHarmonic -""" -def gaussianPacket(x, x0, sigma, p0, extra_param=None): - global hred - return 1./(2*np.pi*sigma**2)**(0.25) * np.exp(-1./4. * ((x-x0)/sigma)**2) * np.exp(1j/hred * p0*(x)) -""" - - -@jit#(cache=True) -def heaviside(x, k=1.): - """ - Heaviside. Analytic approximation: 1/(1 + e^-2kr). Larger k, better approximation, less smooth - """ - return 1./(1. + np.exp(-2.*k*x)) - - -# IMPORTANT -# IMPORTANT -# WARNING: eval() https://nedbatchelder.com/blog/201206/eval_really_is_dangerous.html -# WARNING -# This creates a GLOBAL function. This is why it's allowed to return -# it seems locals can't change, but globals can -# https://stackoverflow.com/questions/41100196/exec-not-working-inside-function-python3-x -from numpy import sin, cos, tan, arcsin, arccos, arctan, hypot, arctan2, degrees, radians, deg2rad, rad2deg,\ - sinh, cosh, tanh, arcsinh, arccosh, arctanh,\ - around, rint, fix, floor, ceil, trunc,\ - exp, expm1, exp2, log, log10, log2, log1p, logaddexp, logaddexp2,\ - power, sqrt, i0, sinc,\ - sign,\ - pi - -ln = log #for commodity - -# Heaviside is not supported in numba - class InvalidFormat(Exception): pass + class MaliciousInput(Exception): pass - def createFunc(expression, variableDict): if expression == "": raise InvalidFormat try: # Things like {px} are substituted by their corresponding actual extra_param expressionFormated = expression.format(**variableDict) except: - raise InvalidFormat #Exception("Could not replace global variables properly") + raise InvalidFormat # Exception("Could not replace global variables properly") ################# # SAFETY CHECKS # ################# - if "print" in expression or "import" in expression or "sys" in expression or "os." in expression or "open" in expression\ - or "__" in expression:# or "__builtins__" in expression: + if "print" in expression or "import" in expression or "sys" in expression or "os." in expression or "open" in expression \ + or "__" in expression: # or "__builtins__" in expression: raise MaliciousInput # This is very very bad. Someone is doing something really wrong - exec(""" def funcManualGLOBAL(x, y, t=0., extra_param=np.array([])): r = sqrt(x*x + y*y) @@ -856,14 +1040,16 @@ def funcManualGLOBAL(x, y, t=0., extra_param=np.array([])): # After defining the function, we test it once, to see if it works. If not, it will raise an exception, # which should be catched where createFunc is used - if not np.isfinite( funcManualGLOBAL(0., 0., 0., np.array([0.5]*100)) ): + if not np.isfinite(funcManualGLOBAL(0., 0., 0., np.array([0.5] * 100))): raise ZeroDivisionError # it can happen for other reasons too though. return funcManualGLOBAL + class WindowManager(ScreenManager): pass + class MainScreen(Screen): def __init__(self, **kwargs): super(Screen, self).__init__(**kwargs) @@ -874,7 +1060,6 @@ class SandboxScreen(Screen): settingsButton = ObjectProperty(None) - def __init__(self, **kwargs): self._first_init() super(SandboxScreen, self).__init__(**kwargs) @@ -885,7 +1070,9 @@ def on_enter(self, *args): self.animation.reset_plot() def _first_init(self): - self.Nx = 200; self.Ny = 200; L = 10. + self.Nx = 200; + self.Ny = 200; + L = 10. self.x0, self.y0 = -L, -L self.xf, self.yf = L, L @@ -899,22 +1086,22 @@ def _first_init(self): self.extra_param[0] = 0. self.extra_param[1] = 5. self.extra_param[2] = 0. - self.extra_param[self.nVar-2] = -5. - self.extra_param[self.nVar-1] = -2.5 + self.extra_param[self.nVar - 2] = -5. + self.extra_param[self.nVar - 1] = -2.5 # A name can be assigned to each of these - self.paramNames = ["Vx"] + ["Vh"] + ["Vw"]+[""]*(self.nVar-5) + ["px"] + ["py"] - self.paramSliders = [True] + [False] + [True] + [False]*(self.nVar-2) + self.paramNames = ["Vx"] + ["Vh"] + ["Vw"] + [""] * (self.nVar - 5) + ["px"] + ["py"] + self.paramSliders = [True] + [False] + [True] + [False] * (self.nVar - 2) self.sliders = [] self.setSliders(firstInit=True) self.variablesDict = {} self.setVarDict() - #self.variablesDict = {'px': 'extra_param[{}]'.format(self.nVar-2), 'py': 'extra_param[{}]'.format(self.nVar-1)} + # self.variablesDict = {'px': 'extra_param[{}]'.format(self.nVar-2), 'py': 'extra_param[{}]'.format(self.nVar-1)} - self.initState = mathPhysics.gaussian2D(7, 1., self.extra_param[self.nVar-2], - 7., 1., self.extra_param[self.nVar-1]) + self.initState = mathPhysics.gaussian2D(7, 1., self.extra_param[self.nVar - 2], + 7., 1., self.extra_param[self.nVar - 1]) self.initStateDef = \ "gaussianPacket(x, 7, 1, {px}) * gaussianPacket(y, 7, 1, {py})" - #"1/(2*pi)**0.5 * exp(-1./4. * ((x-7)**2 + (y-7)**2)) * exp(1j * ({px}*x + {py}*y))" + # "1/(2*pi)**0.5 * exp(-1./4. * ((x-7)**2 + (y-7)**2)) * exp(1j * ({px}*x + {py}*y))" self.potential = mathPhysics.potentialBarrierYCustom self.potentialDef = "{Vh} * exp(-((x-{Vx}) ** 2) / 0.1) / sqrt(0.1 * pi) * (0 if abs(y) < {Vw} else 1)" @@ -925,47 +1112,90 @@ def _first_init(self): self.animation = animate.QuantumAnimation(self.QSystem, dtSim=0.01, dtAnim=0.05, debugTime=True, showPotential=True, varyingPotential=True, - showMomentum=True, showEnergy=True, forceEnergy=True, showNorm=False, forceNorm=True, + showMomentum=True, showEnergy=True, forceEnergy=True, showNorm=False, + forceNorm=True, scalePsi=True, scaleMom=True, isKivy=True, drawClassical=True, - unit_dist=unit_dist, unit_time=unit_time, unit_energy=unit_energy, unit_mom=unit_mom, + unit_dist=unit_dist, unit_time=unit_time, unit_energy=unit_energy, + unit_mom=unit_mom, toolbar=True) self.savedStates = [] self.tempState = {"psi": self.QSystem.psi, "x0": self.QSystem.x0, "xf": self.QSystem.xf - , "y0": self.QSystem.y0, "yf": self.QSystem.yf, - "name": "temp"} + , "y0": self.QSystem.y0, "yf": self.QSystem.yf, + "name": "temp"} def _finish_init(self, dt): self.plotBox = self.ids.plot box = BoxLayout(orientation="vertical") temp = self.animation.fig.canvas._on_size_changed - self.animation.fig.canvas._on_size_changed = lambda *args: None#print("????") - #box.add_widget(self.animation.navigation.actionbar) + self.animation.fig.canvas._on_size_changed = lambda *args: None # print("????") + # box.add_widget(self.animation.navigation.actionbar) nav = self.animation.navigation - gridnav = ColoredGridLayout(cols=9, height=sp(75/2), size_hint_y=None) + gridnav = ColoredGridLayout(cols=11, height=sp(75 / 2), size_hint_y=None) + + mplSize = nav.actionbar.children[0].children[0].size # We make our custom buttons same size as matplotlib's + mplBckCol = (113 / 255., 161 / 255., 179 / 255., + 1) # had to do aproximate it manually #nav.actionbar.children[0].children[0].background_color + mplWhiteCol = (240 / 255, 240 / 255, 240 / 255, 255 / 255) + + # Control how psi is shown. Mod Squared / Phase / Real part / Imaginary part + self.psiDropdown = DropDown() + + options = ["mod2", "phase", "real", "imag"] + + for option in options: + btnPsioption = LightButtonImage(background_normal='', image_src='images/{}.png'.format(option), + background_color=mplWhiteCol, text="", size_hint_y=None, size=mplSize, + size_hint_image=0.8) + btnPsioption.bind( + on_release=lambda btn, option=option: self.animation.reset_plot(psiRepresentation="{}".format(option))) + self.psiDropdown.add_widget(btnPsioption) + + btnPsi = LightButtonImage(background_normal='', image_src='images/psi.png', background_color=mplWhiteCol, + text="", size_hint=(None, 1), size=mplSize, size_hint_image=0.8) + btnPsi.bind(on_release=self.psiDropdown.open) + + gridnav.add_widget(btnPsi) + + btnE = ToggleButton(background_normal='', background_down='', background_color=mplBckCol, state="down", + bold=True, color=(0, 0, 0, 1), text="E", font_size='30sp', size_hint=(None, 1), + size=mplSize) - mplSize = nav.actionbar.children[0].children[0].size # We make our custom buttons same size as matplotlib's - mplBckCol = (113/255., 161/255., 179/255., 1) # had to do aproximate it manually #nav.actionbar.children[0].children[0].background_color - btnE = ToggleButton(background_normal='', background_down='', background_color=mplBckCol, state="down",bold=True, color=(0,0,0,1), text="E", font_size='30sp', size_hint=(None,1), size=mplSize) def showE(*args): self.animation.reset_plot(showEnergy=btnE.state is 'down') - btnE.background_color = mplBckCol if btnE.state is 'down' else (0,0,0,0) + btnE.background_color = mplBckCol if btnE.state is 'down' else (0, 0, 0, 0) + btnE.bind(on_release=showE) gridnav.add_widget(btnE) - btnP = ToggleButton(background_normal='', background_down='', background_color=mplBckCol, state="down", bold=True, color=(0,0,0,1), text="P", font_size='30sp', size_hint=(None,1), size=mplSize) + btnP = ToggleButton(background_normal='', background_down='', background_color=mplBckCol, state="down", + bold=True, color=(0, 0, 0, 1), text="P", font_size='30sp', size_hint=(None, 1), + size=mplSize) + def showP(*args): self.animation.reset_plot(showMomentum=btnP.state is 'down') btnP.background_color = mplBckCol if btnP.state is 'down' else (0, 0, 0, 0) + btnP.bind(on_release=showP) gridnav.add_widget(btnP) + btnV = ToggleButton(background_normal='', background_down='', background_color=mplBckCol, state="down", + bold=True, color=(0, 0, 0, 1), text="V", font_size='30sp', size_hint=(None, 1), + size=mplSize) + + def showV(*args): + self.animation.reset_plot(showPotential=btnV.state is 'down') + btnV.background_color = mplBckCol if btnV.state is 'down' else (0, 0, 0, 0) + + btnV.bind(on_release=showV) + gridnav.add_widget(btnV) + for i in range(9): kid = nav.actionbar.children[0].children[0] nav.actionbar.children[0].remove_widget(kid) - if i != 0 and i != 2: gridnav.add_widget(kid) # save and settings widget don't do anything anyways + if i != 0 and i != 2: gridnav.add_widget(kid) # save and settings widget don't do anything anyways box.add_widget(gridnav) @@ -976,20 +1206,20 @@ def showP(*args): for slider in self.sliders: self.plotBox.add_widget(slider) - def reallowResizing(*args): self.animation.fig.canvas._on_size_changed = temp + Clock.schedule_once(reallowResizing) - #self.ids.renorm.bind(on_release = self.renorm) + # self.ids.renorm.bind(on_release = self.renorm) # Some bug appeared. When initializing everything matplotlib breaks because in generating canvas # it gets a 0 width, and divides by 0. We allow resizing canvas later only, to make sure kivy screens are well defined - #print("yello0") # This would print actually - #Clock.schedule_once(lambda x: print("yello")) # This not. The error happens in between, drawing is scheduled here + # print("yello0") # This would print actually + # Clock.schedule_once(lambda x: print("yello")) # This not. The error happens in between, drawing is scheduled here - #self.settingsButton.bind(on_release = self.dropdown.open) + # self.settingsButton.bind(on_release = self.dropdown.open) global sandbox sandbox = self @@ -1010,8 +1240,18 @@ def renorm(self, dt): print(np.trapz(np.trapz(np.multiply(np.conj(self.QSystem.psi),self.QSystem.psiCopy)))*self.QSystem.dx*self.QSystem.dy)""" self.animation.QSystem.renorm() + """ + The isidentifier() method returns True if the string is a valid identifier, otherwise False. + +A string is considered a valid identifier if it only contains alphanumeric letters (a-z) and (0-9), or underscores (_). A valid identifier cannot start with a number, or contain any spaces. + +""" + def saveState(self): repeated = False + if self.ids.stateName.text == '' or not self.ids.stateName.text.isidentifier(): + TextPopup("Nom invàlid!").open() + return for state in self.savedStates: if self.ids.stateName.text == state["name"]: repeated = True @@ -1020,23 +1260,25 @@ def saveState(self): else: self.savedStates.append({"psi": self.QSystem.psi.copy(), "x0": self.QSystem.x0, "xf": self.QSystem.xf - , "y0": self.QSystem.y0, "yf": self.QSystem.yf, + , "y0": self.QSystem.y0, "yf": self.QSystem.yf, "name": self.ids.stateName.text}) self.ids.stateName.text = "est{}".format(len(self.savedStates)) def setState(self, state): self.QSystem.setState(state) self.animation.manualUpdate(onlyDraw=True) - #self.animation.reset_plot() + # self.animation.reset_plot() def substractComponent(self, state): self.QSystem.substractComponent(state) self.animation.manualUpdate(onlyDraw=True) - #self.animation.reset_plot() + # self.animation.reset_plot() def stopPlaying(self): - try: self.schedule.cancel() - except: pass + try: + self.schedule.cancel() + except: + pass self.ids.pausePlay.state = 'normal' self.paused = True self.animation.paused = True @@ -1052,9 +1294,9 @@ def play(self, dt): self.stopPlaying() else: self.animation.manualUpdate() - #self.animation.update(self.animation.frame) - #self.animation.frame += 1 - #self.animation.fig.canvas.draw() + # self.animation.update(self.animation.frame) + # self.animation.frame += 1 + # self.animation.fig.canvas.draw() def setVarDict(self): self.variablesDict.clear() @@ -1064,9 +1306,9 @@ def setVarDict(self): self.variablesDict = \ {self.paramNames[i]: "extra_param[{}]".format(i) for i in range(self.nVar) if self.paramNames[i] != ""}""" - #print("Variables: ", self.variablesDict) + # print("Variables: ", self.variablesDict) - def setSliders(self, firstInit = False): + def setSliders(self, firstInit=False): if not firstInit: for slider in self.sliders: self.plotBox.remove_widget(slider) @@ -1074,9 +1316,10 @@ def setSliders(self, firstInit = False): self.sliders.clear() for i in range(self.nVar): if self.paramSliders[i]: - self.sliders.append(CustomDataSlider(name=self.paramNames[i], attribute="extra_param", index=i, holder=self, - orientation="vertical", value=float(self.extra_param[i]), - min=float(self.extra_param[i]-10.), max=float(self.extra_param[i]+10.))) + self.sliders.append( + CustomDataSlider(name=self.paramNames[i], attribute="extra_param", index=i, holder=self, + orientation="vertical", value=float(self.extra_param[i]), + min=float(self.extra_param[i] - 10.), max=float(self.extra_param[i] + 10.))) # Careful with numpy, weird interaction with NumericProperty: # https://kivy.org/doc/stable/api-kivy.properties.html#kivy.properties.NumericProperty @@ -1096,15 +1339,13 @@ def newSystem(self): self.animation.resetSystem(self.QSystem) - #self.savedStates.clear() + # self.savedStates.clear() self.tempState = {"psi": self.QSystem.psi, "x0": self.QSystem.x0, "xf": self.QSystem.xf - , "y0": self.QSystem.y0, "yf": self.QSystem.yf, + , "y0": self.QSystem.y0, "yf": self.QSystem.yf, "name": "temp"} - #self.animation.reset_lists() - #self.animation.reset_plot() - - + # self.animation.reset_lists() + # self.animation.reset_plot() class TextPopup(Popup): @@ -1112,58 +1353,64 @@ def __init__(self, text, **kwargs): super(TextPopup, self).__init__(**kwargs) self.add_widget(Label(text=text)) + class PlotPopup(Popup): def __init__(self, data, **kwargs): super(PlotPopup, self).__init__(**kwargs) self.data = data - Clock.schedule_once(self._finish_init) + # self._finish_init() + # Clock.schedule_once(self._finish_init) - def _finish_init(self, *args): - self.plotBox = self.ids.plot + # def _finish_init(self, *args): + plotBox = BoxLayout(padding=20) - self.fig, self.ax = plt.subplots() - FigureCanvasKivyAgg(self.fig) - self.ax = self.fig.gca() # Just in case? - #self.ax = self.fig.add_subplot() # Maybe this doesn't always work + self.fig, ax = plt.subplots() + + # self.ax = self.fig.gca() # Just in case? + # self.ax = self.fig.add_subplot() # Maybe this doesn't always work if self.data["psi"].dtype == np.complex128: - self.psiMod = np.empty((len(self.data["psi"]), len(self.data["psi"][0])), dtype=np.float64) - self.psiMod[:,:] = mathPhysics.abs2(self.data["psi"]) + psiMod = np.empty((len(self.data["psi"]), len(self.data["psi"][0])), dtype=np.float64) + psiMod[:, :] = mathPhysics.abs2(self.data["psi"]) else: - self.psiMod = self.data["psi"] + psiMod = self.data["psi"] # Didn't work on some other computers? ax not defined? - datPlot = self.ax.imshow(self.psiMod.T, origin='lower', - extent=(self.data["x0"], self.data["xf"], self.data["y0"], self.data["yf"]), - aspect = 'equal', cmap = "viridis") + datPlot = ax.imshow(psiMod.T, origin='lower', + extent=(self.data["x0"], self.data["xf"], self.data["y0"], self.data["yf"]), + aspect='equal', cmap="viridis") - self.fig.colorbar(datPlot, ax=self.ax, label=self.data.get("unit_col",r'$(2Å)^{-2}$')) + self.fig.colorbar(datPlot, ax=ax, label=self.data.get("unit_col", r'$(2Å)^{-2}$')) - self.ax.set_xlabel("x ({})".format(self.data.get("unit_ax",'2 Å'))) - self.ax.set_ylabel("y ({})".format(self.data.get("unit_ax", '2 Å'))) + ax.set_xlabel("x ({})".format(self.data.get("unit_ax", '2 Å'))) + ax.set_ylabel("y ({})".format(self.data.get("unit_ax", '2 Å'))) - self.plotBox.add_widget(self.fig.canvas) + plotBox.add_widget(FigureCanvasKivyAgg(self.fig)) + # self.fig.show() + self.add_widget(plotBox) self.fig.canvas.draw() self.fig.canvas.disabled = True - #print(gc.get_referrers(self.fig.canvas)) + # print(gc.get_referrers(self.fig.canvas)) def on_dismiss(self): - super(PlotPopup, self).on_dismiss() - plt.close(self.fig) - # Attempts to fix memory leak, to no avail - self.fig.canvas.disabled = True - self.fig.canvas.canvas.clear() - self.fig.canvas.clear_widgets() - self.plotBox.clear_widgets() - """for atr in vars(self.fig.canvas).copy(): - delattr(self.fig.canvas, atr)""" - delattr(self.fig.canvas, "img_texture") - delattr(self.fig.canvas, "img_rect") - - self.fig.canvas.fig = None - self.fig.canvas = None + # self.fig.canvas.disabled = True + # self.fig.canvas.clear_widgets() + # self.children[0].clear_widgets() + # """for atr in vars(self.fig.canvas).copy(): + # delattr(self.fig.canvas, atr)""" + # delattr(self.fig.canvas, "img_texture") + # delattr(self.fig.canvas, "img_rect") + # + # self.fig.canvas.fig = None + # self.fig.canvas = None + + animate.cleanFigClose(self.fig) + + # del self.fig + # del self.data + super(PlotPopup, self).on_dismiss() class ExamplesScreen(Screen): @@ -1175,33 +1422,93 @@ def __init__(self, **kwargs): def _finish_init(self, *args): # ----------- ----------- ----------- ----------- ----------- # WALL - definition = {"zoomMom":1/3} # Almost default. (None in such case) + definition = {"zoomMom": 1 / 3} # Almost default. (None in such case) self.ids.exampselect.add_widget( Button(text="Barrera", on_release=partial(self.switch, definition=definition))) # ----------- ----------- ----------- ----------- ----------- - # ----------- ----------- ----------- ----------- ----------- # GRAVITY definition = {"initState": mathPhysics.gaussian2D(0., 1., 0., - 6.5, 1., 0.), "potential": mathPhysics.potentialGravity, - "zoomMom":1/3.} + 6.5, 1., 0.), "potential": mathPhysics.potentialGravity, + "zoomMom": 1 / 3., "dtSim": 2 ** (-8), "Ny": 300} self.ids.exampselect.add_widget( Button(text="Gravetat", on_release=partial(self.switch, definition=definition))) # ----------- ----------- ----------- ----------- ----------- + # ----------- ----------- ----------- ----------- ----------- + # FREE (DISPERSION) + + def sliderFree(args, ps): + # box = BoxLayout(orientation='horizontal', size_hint_x=0.3) + return CustomDataSlider(name="σ", attribute="extra_param", index=0, + holder=self.manager.get_screen("playscreen"), + orientation="vertical", min=0.2, max=2, value=None) + # return box + + def initStateFreeDisp(x, y, t=0, extra_param=np.array([2.])): + return gaussianPacket(x, 0., extra_param[0], 0.) * gaussianPacket(y, 0., extra_param[0], 0) + + def freeDispSetup(ax, dat, QSystem, units): + # ax.set_title("Evolució de la dispersió") + ax.set_ylabel(r'$\sigma^2$ ({})'.format(units["unit_dist"])) + ax.set_ylim(top=40) + ax.set_xlabel(r'$t$ ({})'.format(units["unit_time"])) + tlist = np.linspace(0, 15, 300) + sigmaT = QSystem.extra_param[0] ** 2 + tlist ** 2 * ( + mathPhysics.hred / (2 * mathPhysics.M * QSystem.extra_param[0]) ** 2) + dat["teoric"], = ax.plot(tlist, sigmaT, '--') + dat["predict"], = ax.plot(np.array([0.]), np.array([QSystem.varX(0.) ** 2])) + # f = open("checks{:.2f}.dat".format(QSystem.extra_param[0]), "w") + # f.write("") + # f.close() + + def freeDispUpdate(ax, dat, QSystem, units): + lin = dat["predict"] + var = QSystem.extra_param[0] # lin.get_ydata()[0] #!!! This is var^2 !!! + dat["predict"].set_data(np.append(lin.get_xdata(), QSystem.t), + np.append(lin.get_ydata(), QSystem.varX(0.) ** 2)) + # no need to redraw + """# Theoretical checks. This code can be used to compare inner product evolution + x = QSystem.Xmesh + y = QSystem.Ymesh + + psiTeoric = np.sqrt(1./(2.*np.pi*(1j*QSystem.t/2/var + var)**2)) *np.exp(-1./4.*(x**2 + y**2)/(1j/2. * QSystem.t + var**2)) + normaPsi = mathPhysics.euclidNorm(QSystem.psi, QSystem.dx, QSystem.dy) + normaTeo = mathPhysics.euclidNorm(psiTeoric, QSystem.dx, QSystem.dy) + teoricIpsi = mathPhysics.innerProduct2D(psiTeoric, QSystem.psi, QSystem.dx, QSystem.dy)/np.sqrt(normaPsi*normaTeo) + #print(" = ", teoricIpsi/(np.sqrt(normaPsi*normaTeo))) + #print("norma teo : ", normaTeo) + #print("norma psi : ", normaPsi) + #Wow, horrible code + #f = open("checks{:.2f}.dat".format(QSystem.extra_param[0]), "a") + #f.write("{}\t{}\t{}\t{}\n".format(normaPsi,normaTeo,teoricIpsi.real,teoricIpsi.imag)) + #f.close()""" + return False + + definition = {"initState": initStateFreeDisp, "potential": mathPhysics.potential0, + "drawClassical": False, "drawExpected": False, "showEnergy": False, "plotWidget": sliderFree, + 'extra_param': np.array([1.5]), "zoomMom": 0.3, "customPlot": (freeDispSetup, freeDispUpdate), + "step": 'eigen', "dtSim": 0.04 + } + + self.ids.exampselect.add_widget( + Button(text="Dispersió\nLliure", on_release=partial(self.switch, definition=definition))) + # ----------- ----------- ----------- ----------- ----------- # ----------- ----------- ----------- ----------- ----------- # UNCERTAINTY + k = 0.5 + @jit def potentialClosingSoft(x, y, t, extra_param): global L # Heaviside. Analytic approximation: 1/(1 + e^-2kr). Larger k, better approximation r = np.sqrt(x * x + y * y) k = 1 - return 100 * 1 / (1 + np.exp(-2 * k * (r - 10. / 2 + 9.5 / 2 * (1 - 1. / (1 + 0.2 * t))))) # 0.5 originally - + # return 100 * 1 / (1 + np.exp(-2 * k * (r - 10. / 2 + 9.5 / 2 * (1 - 1. / (1 + 0.2 * t))))) # 0.5 originally + return min(100., 1 / 2 * (1 / 4 + 80 * tanh(0.01 * t)) * r ** 2) # +#20*(1-1/(1+0.2*t) ))*r**2) def plotUnc(args, ps): grid = BoxLayout(size_hint_x=0.3, size_hint_y=0.4, padding=0, pos_hint={'center_x': .5, 'center_y': .5}) @@ -1212,10 +1519,11 @@ def setUncertainty(args): args["figUnc"] = plt.figure() args["axUnc"] = plt.axes() args["canvasFig"] = FigureCanvasKivyAggModified(args["figUnc"]) - XUncMin = 1 / np.linspace(0.5,10., 10) # better spacing - YUncMin = mathPhysics.hred/2 / XUncMin + XUncMin = 1 / np.linspace(0.5, 10., 10) # better spacing + YUncMin = mathPhysics.hred / 2 / XUncMin - args["axUnc"].set(xlabel=r'$\sigma_x$ ({})'.format(unit_dist), ylabel=r'$\sigma_{{p_x}}$ ({})'.format(unit_mom), + args["axUnc"].set(xlabel=r'$\sigma_x$ ({})'.format(unit_dist), + ylabel=r'$\sigma_{{p_x}}$ ({})'.format(unit_mom), title='Relació entre incerteses\n(variàncies $\\sigma_{p_x}$ i $\\sigma_x$)') args["axUnc"].plot(XUncMin, YUncMin, 'r--', label=r'$\sigma_x \sigma_{p_x} = \hbar/2$') args["axUnc"].grid() @@ -1224,15 +1532,15 @@ def setUncertainty(args): args["axUnc"].set_xscale('log') args["axUnc"].set_yscale('log') - args["datUnc"], = args["axUnc"].plot(np.array([1.]),np.array([1./2.])) + args["datUnc"], = args["axUnc"].plot(np.array([1.]), np.array([1. / 2.])) args["datUncPoint"], = args["axUnc"].plot(np.array([1.]), np.array([1. / 2.]), 'o') args["sigmax"] = [] - #args["sigmay"] = [] + # args["sigmay"] = [] args["sigmapx"] = [] - #args["sigmapy"] = [] + # args["sigmapy"] = [] def extra_update_unc(args, ps): - #ps = self.manager.get_screen("playscreen") + # ps = self.manager.get_screen("playscreen") varX = ps.QSystem.varX() args["sigmax"].append(varX) varPx = ps.QSystem.varPx() @@ -1241,8 +1549,8 @@ def extra_update_unc(args, ps): args["datUnc"].set_data(args["sigmax"], args["sigmapx"]) args["datUncPoint"].set_data(varX, varPx) - #args["sigmay"].append(ps.QSystem.varY()) - #args["sigmapy"].append(ps.QSystem.varPy()) + # args["sigmay"].append(ps.QSystem.varY()) + # args["sigmapy"].append(ps.QSystem.varPy()) """args["figUnc"].draw_artist(args["axUnc"].patch) args["figUnc"].draw_artist(args["datUnc"]) @@ -1253,47 +1561,86 @@ def extra_on_enter_unc(args): args["figUnc"].tight_layout() args["canvasFig"].draw() + def extra_clean_unc(args): + animate.cleanFigClose(args["figUnc"]) + def uncertaintyPlotSetup(ax, dat, QSystem, units, zoomOut=False): + XUncMin = 1 / np.linspace(0.5, 10., 10) # better spacing + YUncMin = mathPhysics.hred / 2 / XUncMin + ax.set(xlabel=r'$\sigma_x$ ({})'.format(units["unit_dist"]), + ylabel=r'$\sigma_{{p_x}}$ ({})'.format(units["unit_mom"]), + title='Relació entre incerteses\n(variàncies $\\sigma_{p_x}$ i $\\sigma_x$)') - def extra_clean_unc(args): - plt.close(args["figUnc"]) + zoom = 3 if zoomOut else 1 + ax.set_xlim(XUncMin[-1], XUncMin[0] * zoom) + ax.set_ylim(YUncMin[0], YUncMin[-1] * zoom) + + ax.plot(XUncMin, YUncMin, 'r--', label=r'$\sigma_x \sigma_{p_x} = \hbar/2$') + + # ax.grid() # Problems with grid disappearing when redrawing plot (do delete trace of points) + + ax.legend() + + ax.set_xscale('log') + ax.set_yscale('log') + + dat["datUnc"], = ax.plot(np.array([1.]), np.array([1. / 2.])) + dat["datUncPoint"], = ax.plot(np.array([1.]), np.array([1. / 2.]), 'o') + dat["sigmax"] = [] + dat["sigmapx"] = [] + + def uncertaintyPlotUpdate(ax, dat, QSystem, units): + varX = QSystem.varX() + dat["sigmax"].append(varX) + varPx = QSystem.varPx() + dat["sigmapx"].append(varPx) + + """for it, patch in enumerate(ax.patches): + dat["p{}".format(it)] = patch""" + + dat["datUnc"].set_data(dat["sigmax"], dat["sigmapx"]) + dat["datUncPoint"].set_data(varX, varPx) + + ax.redraw_in_frame() definition = {"initState": mathPhysics.gaussian2D(0., 1., 0., 0., 1., 0.), "potential": potentialClosingSoft, - "drawClassical":False, "drawExpected": False, "showEnergy":False, "plotWidget": plotUnc, "zoomMom":0.6, - "extra_update": extra_update_unc, "extra_clean": extra_clean_unc, "extra_on_enter": extra_on_enter_unc - } + "drawClassical": False, "drawExpected": False, "showEnergy": False, "zoomMom": 0.6, + "customPlot": (uncertaintyPlotSetup, uncertaintyPlotUpdate), "customPlotUpdate": True} + # "extra_update": extra_update_unc, "extra_clean": extra_clean_unc, "extra_on_enter": extra_on_enter_unc, "plotWidget": plotUnc, + # } doubleUncertainty = GridLayout(rows=2) uncAutom = Button(text="Principi d'incertesa", on_release=partial(self.switch, definition=definition, setExtraArgs=setUncertainty)) - @jit def potentialClosingManual(x, y, t, extra_param): global L # Heaviside. Analytic approximation: 1/(1 + e^-2kr). Larger k, better approximation r = np.sqrt(x * x + y * y) - return 100 * heaviside(r-extra_param[0], 1.) + return 100 * heaviside(r - extra_param[0], 1.) def plotUncManual(args, ps): - box = BoxLayout(orientation='horizontal', size_hint_x=0.3) - box.add_widget(plotUnc(args, ps)) - box.add_widget(CustomDataSlider(name="R", attribute="extra_param", index=0, holder=self.manager.get_screen("playscreen"), - orientation="vertical", min=0., max=10., value=5.)) - return box + # box = BoxLayout(orientation='horizontal')#, size_hint_x=0.3) + # box.add_widget(plotUnc(args, ps)) + # box.add_widget( + return CustomDataSlider(name="R", attribute="extra_param", index=0, + holder=self.manager.get_screen("playscreen"), + orientation="vertical", min=0., max=10., value=5.) + # return box definition = {"initState": mathPhysics.gaussian2D(0., 1., 0., 0., 1., 0.), "potential": potentialClosingManual, - "drawClassical": False, "drawExpected": False, "showEnergy": False, "plotWidget": plotUncManual, - "extra_update": extra_update_unc, "extra_clean": extra_clean_unc, - "extra_on_enter": extra_on_enter_unc, 'extra_param':np.array([5.]), "zoomMom":0.6 - } - + "drawClassical": False, "drawExpected": False, "showEnergy": False, 'extra_param': np.array([5.]), + "zoomMom": 0.6, + "customPlot": (partial(uncertaintyPlotSetup, zoomOut=True), uncertaintyPlotUpdate), + "customPlotUpdate": True, "plotWidget": plotUncManual + } # "extra_update": extra_update_unc, "extra_clean": extra_clean_unc, "extra_on_enter": extra_on_enter_unc, uncManual = Button(text="Principi d'incertesa (manual)", - on_release=partial(self.switch, definition=definition, setExtraArgs=setUncertainty)) + on_release=partial(self.switch, definition=definition, setExtraArgs=setUncertainty)) doubleUncertainty.add_widget(uncAutom) doubleUncertainty.add_widget(uncManual) @@ -1306,9 +1653,12 @@ def plotUncManual(args, ps): 3., 0.6, 0.), "potential": mathPhysics.potentialHarmonic, "extra_param": np.array([2.]), - "dtSim": 0.005, "scalePot":False, "zoomMom":1/3., "Nx": 300, "Ny": 300, # dtSim: 0.00390625 2**(-8) - "drawClassical":True, "drawClassicalTrace":True, "drawExpected":True, "drawExpectedTrace":True, - "plotWidget": CustomDataSlider(name="k", attribute="extra_param", index=0, holder=self.manager.get_screen("playscreen"), + "dtSim": 0.005, "scalePot": False, "zoomMom": 1 / 3., "Nx": 300, "Ny": 300, + # dtSim: 0.00390625 2**(-8) + "drawClassical": True, "drawClassicalTrace": True, "drawExpected": True, + "drawExpectedTrace": True, + "plotWidget": CustomDataSlider(name="k", attribute="extra_param", index=0, + holder=self.manager.get_screen("playscreen"), orientation="vertical", min=1., max=2., value=2.)} self.ids.exampselect.add_widget( Button(text="Oscil·lador Harmònic", on_release=partial(self.switch, definition=definition))) @@ -1318,24 +1668,27 @@ def plotUncManual(args, ps): # ----------- ----------- ----------- ----------- ----------- # Double Slit - @jit def slit(x, n, width, dist): - return x <= n/2*width + dist*(n-1)/2 and (x-width/2*(n%2)+dist/2*((n+1)%2))/(dist+width) % 1 >= dist/(dist+width) + # only defined for positive x + return x <= n / 2 * width + dist * (n - 1) / 2 and (x - width / 2 * (n % 2) + dist / 2 * ((n + 1) % 2)) / ( + dist + width) % 1 >= dist / (dist + width) - #----[ ]--[ ]--[ ]---- - #------[ ]---[ ]------ + # ----[ ]--[ ]--[ ]---- + # ------[ ]---[ ]------ @jit def potentialDoubleSlit(x, y, t, extra_param): # extra_param: [nSlits, slitWidth, slitSeparation, wallWidth] - return 400. if (abs(x) 10/p0A: - instance.datPsi.set_clim(vmax=np.max(instance.QSystem.psiMod[int(10./instance.QSystem.dx):].T), vmin=0.) + def extraUpdateClimSlit(instance): # We want to clearly see what propagates + if instance.QSystem.t > 10 / p0A: + instance.datPsi.set_clim(vmax=np.max(instance.QSystem.psiMod[int(10. / instance.QSystem.dx):].T), + vmin=0.) ps = self.manager.get_screen("playscreen") Qs = ps.QSystem @@ -1413,30 +1774,38 @@ def extraUpdateClimSlit(instance): # We want to cle definition = {"initState": mathPhysics.gaussian2D(-4., 1, p0A, 0., 2, 0.), - "potential": potentialDoubleSlit, "extraUpdates": [extraUpdateClimSlit], "extraUpdatesStart": True, - "extra_param": np.array([2, 1, 1, 0.5, 0., 0.]), "x0": x0A, "xf": xfA, "Nx": NxA, "Ny":NxA, "y0":-10, "yf": 10, - "dtSim": 2**(-7), "stepsPerFrame":2, "scalePot": False, + "potential": potentialDoubleSlit, "extraUpdates": [extraUpdateClimSlit], + "extraUpdatesStart": True, + "extra_param": np.array([2, 1, 1, 0.5, 0., 0.]), "x0": x0A, "xf": xfA, "Nx": NxA, "Ny": NyA, + "y0": -10, "yf": 10, + "dtSim": dtSimA, "stepsPerFrame": 2, "scalePot": False, "drawClassical": False, "drawExpected": False, "customPlot": (slitInterferenceSetup, slitInterferenceUpdate), "customPlotFull": True, - "showMomentum": False, "showEnergy": False, "showNorm": False, "duration": 1., "plotWidget": slidersSlit} + "showMomentum": False, "showEnergy": False, "showNorm": False, "duration": 1., + "plotWidget": slidersSlit} self.ids.exampselect.add_widget( Button(text="Doble Escletxa", on_release=partial(self.switch, definition=definition))) # ----------- ----------- ----------- ----------- ----------- # Double Slit + Aharonov Bohm effect # Very similar to double slit, just some extra things + # MEANING OF ALPHA = -e FLUX / hc, which is numerically 2pi e FLUX/hred c. Which is [2pi FLUX]!!! phase (real: e/hred · flux ) def slidersAharonov(*args): - box = BoxLayout(orientation='horizontal', width=sp(80/2)*3, size_hint_x=None) - stack = BoxLayout(orientation='vertical', width=sp(80/2), spacing=sp(10)) - stack.add_widget(CustomDataSlider(name="n", attribute="extra_param", index=0, holder=self.manager.get_screen("playscreen"), - orientation="vertical", min=2, max=10, step=2, value=None, variableLimits=False)) + box = BoxLayout(orientation='horizontal', width=sp(80 / 2) * 3, size_hint_x=None) + stack = BoxLayout(orientation='vertical', width=sp(80 / 2), spacing=sp(10)) + stack.add_widget(CustomDataSlider(name="n", attribute="extra_param", index=0, + holder=self.manager.get_screen("playscreen"), + orientation="vertical", min=2, max=10, step=2, value=None, + variableLimits=False)) stack.add_widget(CustomDataSlider(name="w", attribute="extra_param", index=1, - holder=self.manager.get_screen("playscreen"), - orientation="vertical", min=0., max=3., step=0.1, value=None, variableLimits=False)) + holder=self.manager.get_screen("playscreen"), + orientation="vertical", min=0., max=3., step=0.1, value=None, + variableLimits=False)) stack.add_widget(CustomDataSlider(name="d", attribute="extra_param", index=2, - holder=self.manager.get_screen("playscreen"), - orientation="vertical", min=0.5, max=3., step=0.1, value=None, variableLimits=False)) + holder=self.manager.get_screen("playscreen"), + orientation="vertical", min=0.5, max=3., step=0.1, value=None, + variableLimits=False)) box.add_widget(stack) box.add_widget(CustomDataSlider(name="α", attribute="extra_param", index=4, holder=self.manager.get_screen("playscreen"), @@ -1447,13 +1816,11 @@ def slidersAharonov(*args): return box - - def AharonovA(x, y, alpha): - k = alpha/(x*x+y*y+1e-10) + x = x - 1 + k = alpha / (x * x + y * y + 1e-10) return k * -y, k * x - def aharonovExtraDrawings(instance=None): ps = self.manager.get_screen("playscreen") Qs = ps.QSystem @@ -1466,19 +1833,24 @@ def aharonovExtraDrawings(instance=None): if instance.firstDraw: ps.extraArgs["vectorPotential"] = instance.axPsi.quiver(Xm.T, Ym.T, Amesh[0].T, Amesh[1].T, pivot="mid", color='red', alpha=0.2, scale=5) - elif instance.updating: ps.extraArgs["vectorPotential"].set_UVC(Amesh[0].T, Amesh[1].T) + elif instance.updating: + ps.extraArgs["vectorPotential"].set_UVC(Amesh[0].T, Amesh[1].T) return ps.extraArgs["vectorPotential"] # A ~ alpha/r^2 (-y, x, 0) definition = {"initState": mathPhysics.gaussian2D(-4., 1, p0A, 0., 2, 0.), - "potential": potentialDoubleSlit, "extraUpdates": [extraUpdateClimSlit, aharonovExtraDrawings], "extraUpdatesStart":True, "extraUpdatesUpdate":True, - "extra_param": np.array([2, 1, 1, 0.5, 0.4, 0]), "x0": x0A, "xf": xfA, "Nx": NxA, "Ny": NxA, "y0": -10, "yf": 10, - "dtSim": 2 ** (-7), "stepsPerFrame":2, "scalePot": False, + "potential": potentialDoubleSlit, "extraUpdates": [extraUpdateClimSlit, aharonovExtraDrawings], + "extraUpdatesStart": True, "extraUpdatesUpdate": True, + "extra_param": np.array([2, 1, 1, 0.5, 0.4, 0]), "x0": x0A, "xf": xfA, "Nx": NxA, "Ny": NyA, + "y0": -10, "yf": 10, + "dtSim": dtSimA, "stepsPerFrame": 2, "scalePot": False, "drawClassical": False, "drawExpected": False, - "customPlot": (slitInterferenceSetup, slitInterferenceUpdate), "customPlotUpdate":True, "customPlotFull":True, + "customPlot": (slitInterferenceSetup, slitInterferenceUpdate), "customPlotUpdate": True, + "customPlotFull": True, "showMomentum": False, "showEnergy": False, "showNorm": False, "duration": 1., + "psiRepresentation": "phase", "plotWidget": slidersAharonov, "customOperator": mathPhysics.aharonovBohmOperator, } @@ -1492,13 +1864,9 @@ def switch(self, *args, definition=None, setExtraArgs=None): self.manager.get_screen("playscreen").sourceScreen = "examples" - - - - class GameCoin(Widget): - def __init__(self, callbackClick = None, callbackMiss = None, callbackEnd = None, duration=2.5, **kwargs): - super(GameCoin, self).__init__(width=dp(50/2), height=dp(50/2), **kwargs) + def __init__(self, callbackClick=None, callbackMiss=None, callbackEnd=None, duration=2.5, **kwargs): + super(GameCoin, self).__init__(width=dp(50 / 2), height=dp(50 / 2), **kwargs) self.cbClick = callbackClick self.cbMiss = callbackMiss self.cbEnd = callbackEnd @@ -1513,6 +1881,7 @@ def complete(*args): if self.parent is not None: self.parent.remove_widget(self) if self.cbMiss is not None: self.cbMiss() if self.cbEnd is not None: self.cbEnd() + self.anim = Animation(opacity=0, duration=self.duration, t="in_quint") self.anim.bind(on_complete=complete) self.anim.start(self.ids.image) @@ -1527,8 +1896,9 @@ def on_touch_down(self, touch): def complete(*args): if self.parent is not None: self.parent.remove_widget(self) if self.cbEnd is not None: self.cbEnd() - #little jump - self.anim = Animation(pos_hint={"center_x":1.2}, duration=0.6)# + Animation(duration=0.2) + + # little jump + self.anim = Animation(pos_hint={"center_x": 1.2}, duration=0.6) # + Animation(duration=0.2) self.anim.bind(on_complete=complete) self.anim.start(self) return super(GameCoin, self).on_touch_down(touch) @@ -1549,8 +1919,8 @@ def _finish_init(self, *args): @njit([numba.float64(numba.float64, numba.float64, numba.float64, numba.float64[:])]) def potentialHarmonicWellMovingSoft(x, y, t, extra_param): res = 1 / 2. * extra_param[5] * ((x - extra_param[0] - extra_param[2] * (t - extra_param[4])) ** 2 - + (y - extra_param[1] - extra_param[3] * (t - extra_param[4])) ** 2 - ) + + (y - extra_param[1] - extra_param[3] * (t - extra_param[4])) ** 2 + ) if res > extra_param[6]: return extra_param[6] return res @@ -1559,8 +1929,8 @@ def potentialHarmonicWellMoving(x, y, t, extra_param): res = 1 / 2. * extra_param[5] * ((x - extra_param[2]) ** 2 + (y - extra_param[3]) ** 2) if res > extra_param[6]: return extra_param[6] return res - inicial2D = mathPhysics.eigenvectorHarmonic2DGenerator(0., 1, 0., 1, kHarm) # 2 2 + inicial2D = mathPhysics.eigenvectorHarmonic2DGenerator(0., 0, 0., 0, kHarm) # 2 2 # Leaderboard stuff (json file) # lazily copied from: https://stackoverflow.com/a/12309296 def json_load_file(json_file): @@ -1575,7 +1945,7 @@ def json_dump_to_file(json_file, json_dict): class LeaderboardInput(TextInput): # Only 5 characters, and upper leter def insert_text(self, substring, from_undo=False): - s = substring.upper()[:5-len(self.text)] + s = substring.upper()[:5 - len(self.text)] return super().insert_text(s, from_undo=from_undo) def setMoveGame(args): @@ -1587,10 +1957,13 @@ def setMoveGame(args): args["score"] = 0 args["health"] = 3 args["coins"] = 0 - args["newCoin"]=True # New coin can be scheduled to be created - args["goalX"] = ps.QSystem.x0 + args["goalRadius"] + random.random() * (ps.QSystem.xf - ps.QSystem.x0 - 2 * args["goalRadius"]) - args["goalY"] = ps.QSystem.y0 + args["goalRadius"] + random.random() * (ps.QSystem.yf - ps.QSystem.y0 - 2 * args["goalRadius"]) - args["goalCircle"] = plt.Circle((args["goalX"], args["goalY"]), args["goalRadius"], alpha=0.2, color='black') + args["newCoin"] = True # New coin can be scheduled to be created + args["goalX"] = ps.QSystem.x0 + args["goalRadius"] + random.random() * ( + ps.QSystem.xf - ps.QSystem.x0 - 2 * args["goalRadius"]) + args["goalY"] = ps.QSystem.y0 + args["goalRadius"] + random.random() * ( + ps.QSystem.yf - ps.QSystem.y0 - 2 * args["goalRadius"]) + args["goalCircle"] = plt.Circle((args["goalX"], args["goalY"]), args["goalRadius"], alpha=0.2, + color='black') args["firstDraw"] = True args["drawnCircle"] = None args["gameOver"] = False @@ -1604,7 +1977,6 @@ def setMoveGame(args): # - With keyboard # - Directly with mouse!? - def extra_keyboard_movegame(ps, keyboard, keycode, text, modifiers): if keycode[1] == 'up' or 'down' or 'left' or 'right': t = ps.QSystem.t @@ -1621,18 +1993,18 @@ def extra_keyboard_movegame(ps, keyboard, keycode, text, modifiers): ps.extra_param[2] += 0.25 return True - - def drawCircle(instance=None): #Game Logic. Done in matplotlib directly (to be able to draw circles there), - # which makes it a bit more of a mess than it needs to be, but it's the same + def drawCircle(instance=None): # Game Logic. Done in matplotlib directly (to be able to draw circles there), + # which makes it a bit more of a mess than it needs to be, but it's the same ps = self.manager.get_screen("playscreen") - ps.extraArgs["progress"].value = (instance.frame%(6/0.04)) / (6/0.04) * 100 + ps.extraArgs["progress"].value = (instance.frame % (6 / 0.04)) / (6 / 0.04) * 100 if ps.extraArgs["energyButton"].disabled: fig = ps.extraArgs["energyGraph"] # Redraw background. If not, bars don't get deleted when enregy goes down for patch in ps.extraArgs["energyAx"].patches: - if patch != ps.extraArgs["energyK"][0] and patch != ps.extraArgs["energyV"][0]: fig.draw_artist(patch) + if patch != ps.extraArgs["energyK"][0] and patch != ps.extraArgs["energyV"][0]: fig.draw_artist( + patch) EKinetic = np.real(ps.QSystem.kineticEnergy()) for bar in ps.extraArgs["energyK"]: @@ -1649,47 +2021,50 @@ def drawCircle(instance=None): #Game Logic. Done in matplotlib directly (to be instance.drawnCircle = instance.axPsi.add_patch(ps.extraArgs["goalCircle"]) ps.extraArgs["firstDraw"] = False instance.observedParticle = None - ps.playButton.disabled = True # Pausing can cause cheating. L - elif instance.frame%(6/0.04)==0: # Every 6 seconds + ps.playButton.disabled = True # Pausing can cause cheating. L + elif instance.frame % (6 / 0.04) == 0: # Every 6 seconds instance.QSystem.modSquared() i, j = mathPhysics.generateAsDiscreteDistribution(instance.QSystem.psiMod) x, y = instance.QSystem.X[i], instance.QSystem.X[j] instance.observedParticle = instance.axPsi.add_patch(plt.Circle((x, y), radius=0.15, color='cyan')) - if (x-ps.extraArgs["goalX"])**2 + (y-ps.extraArgs["goalY"])**2 <= ps.extraArgs["goalRadius"]**2: + if (x - ps.extraArgs["goalX"]) ** 2 + (y - ps.extraArgs["goalY"]) ** 2 <= ps.extraArgs[ + "goalRadius"] ** 2: instance.drawnCircle.set(color='green', alpha=0.5) - ps.extraArgs["score"] += 5 * ps.extraArgs["difficulty"] + 2 * (ps.extraArgs["difficulty"] - 1)**2 + ps.extraArgs["score"] += 5 * ps.extraArgs["difficulty"] + 2 * (ps.extraArgs["difficulty"] - 1) ** 2 ps.extraArgs["labelScore"].text = "Punts: {}".format(ps.extraArgs["score"]) else: instance.drawnCircle.set(color='red', alpha=0.5) ps.extraArgs["health"] -= 1 - ps.extraArgs["labelHealth"].text="Vides: "+'♥' * ps.extraArgs["health"] + ps.extraArgs["labelHealth"].text = "Vides: " + '♥' * ps.extraArgs["health"] if ps.extraArgs["health"] == 0: ps.stopPlaying() - #ps.playButton.disabled_color = 'red' + # ps.playButton.disabled_color = 'red' ps.playButton.disabled = True ps.extraArgs["gameOver"] = True ps.extraArgs["shop"].disabled = True lb = ps.extraArgs["leaderboard"] - ranking = 6 # By default, not winner - for pos in range(1,5+1): + ranking = 6 # By default, not winner + for pos in range(1, 5 + 1): if lb[str(pos)]["score"] < ps.extraArgs["score"]: ranking = pos break if ranking != 6: for pos in range(5, ranking, -1): - lb[str(pos)] = lb[str(pos-1)] - - lb[str(ranking)] = {"name":'', "score":ps.extraArgs["score"]} + lb[str(pos)] = lb[str(pos - 1)] + lb[str(ranking)] = {"name": '', "score": ps.extraArgs["score"]} layout = BoxLayout(orientation='vertical') layout.add_widget(Label(text="Bona puntuació!\nIdentifica't:")) - nameInput = LeaderboardInput(multiline=False, font_name='RobotoMono-Regular', size_hint=(None,None), height=50, width=100, pos_hint={"center_x":0.5, "center_y":0.5}) + nameInput = LeaderboardInput(multiline=False, font_name='RobotoMono-Regular', + size_hint=(None, None), height=50, width=100, + pos_hint={"center_x": 0.5, "center_y": 0.5}) + def nameSet(instance): lb[str(ranking)]["name"] = nameInput.text json_dump_to_file(ps.extraArgs["lbDir"], lb) @@ -1697,41 +2072,46 @@ def nameSet(instance): lbpopup = Popup(size_hint=(0.4, 0.4), title='Ranking') lbpopup.add_widget(Label( text="{:>5}{:>10}{:>10}\n{}\n".format("POS", "NOM", "PUNTS", "-" * 25) + "".join( - ["{:>5}{:>10}{:>10}\n".format(num, lb[str(num)]["name"], lb[str(num)]["score"]) for num in + ["{:>5}{:>10}{:>10}\n".format(num, lb[str(num)]["name"], lb[str(num)]["score"]) + for num in lb]), font_name='RobotoMono-Regular')) # Font comes with kivy, so for sure everyone has it lbpopup.open() layout.add_widget(nameInput) - newRecord = Popup(size_hint=(0.4,0.4), title='Record!', on_dismiss=nameSet) + newRecord = Popup(size_hint=(0.4, 0.4), title='Record!', on_dismiss=nameSet) newRecord.add_widget(layout) newRecord.open() else: - lbpopup = Popup(size_hint=(0.4,0.4), title='Ranking') - lbpopup.add_widget(Label(text="{:>5}{:>10}{:>10}\n{}\n".format("POS", "NOM", "PUNTS", "-" * 25) + "".join( - ["{:>5}{:>10}{:>10}\n".format(num, lb[str(num)]["name"], lb[str(num)]["score"]) for num in lb]),font_name='RobotoMono-Regular')) # Font comes with kivy, so for sure everyone has it + lbpopup = Popup(size_hint=(0.4, 0.4), title='Ranking') + lbpopup.add_widget( + Label(text="{:>5}{:>10}{:>10}\n{}\n".format("POS", "NOM", "PUNTS", "-" * 25) + "".join( + ["{:>5}{:>10}{:>10}\n".format(num, lb[str(num)]["name"], lb[str(num)]["score"]) for + num in lb]), + font_name='RobotoMono-Regular')) # Font comes with kivy, so for sure everyone has it lbpopup.open() return instance.axPsi.text(0.5, 0.5, 'GAME OVER!', dict(size=30, fontweight=800, color='white'), - horizontalalignment='center', verticalalignment='center', - path_effects=[animate.peffects.withStroke(linewidth=4, foreground="black")], - transform=instance.axPsi.transAxes), instance.drawnCircle, instance.observedParticle + horizontalalignment='center', verticalalignment='center', + path_effects=[ + animate.peffects.withStroke(linewidth=4, foreground="black")], + transform=instance.axPsi.transAxes), instance.drawnCircle, instance.observedParticle ps.stopPlaying() def newStep(*args): ps.extraArgs["goalRadius"] = ps.extraArgs["goalRadius"] ** 0.9 ps.extraArgs["goalX"] = instance.QSystem.x0 + ps.extraArgs["goalRadius"] + random.random() * ( - instance.QSystem.xf - instance.QSystem.x0 - 2 * ps.extraArgs["goalRadius"]) + instance.QSystem.xf - instance.QSystem.x0 - 2 * ps.extraArgs["goalRadius"]) ps.extraArgs["goalY"] = instance.QSystem.y0 + ps.extraArgs["goalRadius"] + random.random() * ( - instance.QSystem.yf - instance.QSystem.y0 - 2 * ps.extraArgs["goalRadius"]) + instance.QSystem.yf - instance.QSystem.y0 - 2 * ps.extraArgs["goalRadius"]) instance.drawnCircle.remove() ps.extraArgs["goalCircle"] = plt.Circle((ps.extraArgs["goalX"], ps.extraArgs["goalY"]), ps.extraArgs["goalRadius"], alpha=0.2, color='black') instance.drawnCircle = instance.axPsi.add_patch(ps.extraArgs["goalCircle"]) instance.observedParticle.remove() - #print(ps.extraArgs["score"]) + # print(ps.extraArgs["score"]) ps.startPlaying() Clock.schedule_once(newStep, timeout=1) @@ -1746,7 +2126,7 @@ def newStep(*args): """def extra_update_move(args):""" def coins_add(args, val): - args["coins"]+=val + args["coins"] += val args["labelCoins"].text = "Monedes: {}".format(args["coins"]) # dt is irrelevant, but is an argument used by schedule_once @@ -1759,49 +2139,59 @@ def coin_end(*arg): args["newCoin"] = True ps.plotRelative.add_widget(GameCoin(callbackClick=coin_get, callbackEnd=coin_end, - pos_hint={"center_x":0.2+random.random()*0.6, "center_y":0.2+random.random()*0.6})) - + pos_hint={"center_x": 0.2 + random.random() * 0.6, + "center_y": 0.2 + random.random() * 0.6})) def extra_update_movement(args, ps): if args["newCoin"] is True: args["newCoin"] = False - args["coinSchedule"] = Clock.schedule_once(lambda dt: partial(createCoin, args, ps)(), 5/args["difficulty"]) # Too fast??? Clock.schedule_once doesn't work too well + args["coinSchedule"] = Clock.schedule_once(lambda dt: partial(createCoin, args, ps)(), 5 / args[ + "difficulty"]) # Too fast??? Clock.schedule_once doesn't work too well def extra_info_movement(args, ps): - layout = GridLayout(rows=4,cols=1, width = sp(500/2), size_hint_x=None, padding=20) - - args["labelScore"] = Label(text="Punts: {}".format(args["score"]), size_hint_min_y=sp(22), size_hint=(1,0.07)) - try: args["labelHealth"] = Label(text="Vides: " + '♥'*args["health"], font_name='Arial', color=(1,0,0,1), size_hint_min_y=sp(22), size_hint=(1,0.07)) #Maybe Arial not in all machines? - except: args["labelHealth"] = Label(text="Vides: {}".format(args["health"]), size_hint_min_y=sp(30), size_hint=(1,0.07)) # No heart emoji then + layout = GridLayout(rows=4, cols=1, width=sp(500 / 2), size_hint_x=None, padding=20) + + args["labelScore"] = Label(text="Punts: {}".format(args["score"]), size_hint_min_y=sp(22), + size_hint=(1, 0.07)) + try: + args["labelHealth"] = Label(text="Vides: " + '♥' * args["health"], font_name='Arial', + color=(1, 0, 0, 1), size_hint_min_y=sp(22), + size_hint=(1, 0.07)) # Maybe Arial not in all machines? + except: + args["labelHealth"] = Label(text="Vides: {}".format(args["health"]), size_hint_min_y=sp(30), + size_hint=(1, 0.07)) # No heart emoji then layout.add_widget(args["labelScore"]) layout.add_widget(args["labelHealth"]) - args["labelCoins"] = Label(text="Monedes: {}".format(args["coins"]), color=(0.7,0.6,0,1), size_hint_min_y=sp(22), size_hint=(1,0.07)) + args["labelCoins"] = Label(text="Monedes: {}".format(args["coins"]), color=(0.7, 0.6, 0, 1), + size_hint_min_y=sp(22), size_hint=(1, 0.07)) layout.add_widget(args["labelCoins"]) args["coinSchedule"] = None - args["progress"] = ProgressBar(size_hint_y=None, height=dp(25/2)) - #layout.add_widget(args["progress"]) + args["progress"] = ProgressBar(size_hint_y=None, height=dp(25 / 2)) + # layout.add_widget(args["progress"]) ps.plotRelative.add_widget(args["progress"]) BtnHeight = sp(40) shopBig = BoxLayout(orientation="vertical") args["shop"] = shopBig - shopBig.add_widget(Label(text="________________________\n--------------- Tenda ---------------", size_hint_min_y=sp(60), size_hint=(1,0.14))) - shop = GridLayout(rows=2, cols=2, spacing=sp(5))#, size_hint_y=0.3) + shopBig.add_widget( + Label(text="________________________\n--------------- Tenda ---------------", size_hint_min_y=sp(60), + size_hint=(1, 0.14))) + shop = GridLayout(rows=2, cols=2, spacing=sp(5)) # , size_hint_y=0.3) def buyEnergy(btn): if args["coins"] >= 1: coins_add(args, -1) - args["energyButton"].disabled = True - energyGrid = GridLayout(cols=2,rows=1)#,size_hint_y=0.6) - energyButtonBox = BoxLayout(padding=(sp(10),sp(20))) + energyGrid = GridLayout(cols=2, rows=1) # ,size_hint_y=0.6) + energyButtonBox = BoxLayout(padding=(sp(10), sp(20))) buttonPad = BoxLayout() - args["energyButton"] = Button(text = "(1)Mostra\nEnergia", on_release=buyEnergy, size_hint_y=None, height=BtnHeight) + args["energyButton"] = Button(text="(1)Mostra\nEnergia", on_release=buyEnergy, size_hint_y=None, + height=BtnHeight) buttonPad.add_widget(args["energyButton"]) energyButtonBox.add_widget(buttonPad) @@ -1811,11 +2201,11 @@ def buyEnergy(btn): fig = args["energyGraph"] ax = fig.gca() args["energyAx"] = ax - #ax.set_xlabel("") + # ax.set_xlabel("") ax.get_xaxis().set_visible(False) ax.set_ylabel("Energia ({})".format(unit_energy)) ax.set_ylim([0., 100.]) - #ax.set_xticklabels([]) + # ax.set_xticklabels([]) plt.tick_params( axis='x', # changes apply to the x-axis which='both', # both major and minor ticks are affected @@ -1824,13 +2214,13 @@ def buyEnergy(btn): labelbottom=False) # labels along the bottom edge are off ax.spines['top'].set_visible(False) ax.spines['right'].set_visible(False) - args["energyK"] = ax.bar([0.],[np.real(ps.QSystem.kineticEnergy())], color='b') + args["energyK"] = ax.bar([0.], [np.real(ps.QSystem.kineticEnergy())], color='b') args["energyV"] = ax.bar([0.], [np.real(ps.QSystem.potentialEnergy())], bottom=[1.], color='y') - ax.bar([0.], [500.], color='black') # "Background" - #ax.legend(["K","V"], bbox_to_anchor=(1.04,1)) + ax.bar([0.], [500.], color='black') # "Background" + # ax.legend(["K","V"], bbox_to_anchor=(1.04,1)) FigureCanvasKivyAggModified(fig) - #fig.canvas.size_hint_y=None - #fig.canvas.height=300 + # fig.canvas.size_hint_y=None + # fig.canvas.height=300 energyGrid.add_widget(fig.canvas) shopBig.add_widget(energyGrid) @@ -1838,20 +2228,22 @@ def buyEnergy(btn): def prep(dt): fig.tight_layout() fig.canvas.draw() + Clock.schedule_once(prep) - #Buy health, button + # Buy health, button def extra_life(btn): if args["coins"] >= args["healthCost"]: coins_add(args, -args["healthCost"]) - args["healthCost"]+=2 + args["healthCost"] += 2 args["healthButton"].text = "({})Vida\nExtra".format(args["healthCost"]) args["health"] += 1 args["labelHealth"].text = "Vides: " + '♥' * args["health"] args["healthCost"] = 2 - args["healthButton"] = Button(text="({})Vida\nExtra".format(args["healthCost"]), on_release=extra_life, size_hint_y=None, height=BtnHeight) + args["healthButton"] = Button(text="({})Vida\nExtra".format(args["healthCost"]), on_release=extra_life, + size_hint_y=None, height=BtnHeight) shop.add_widget(args["healthButton"]) # Convert coins to points: @@ -1861,48 +2253,55 @@ def coin_to_point(btn): args["score"] += 1 args["labelScore"].text = "Punts: {}".format(ps.extraArgs["score"]) - shop.add_widget(Button(text="(2)Compra\nPunt", on_release=coin_to_point, size_hint_y=None, height=BtnHeight)) + shop.add_widget( + Button(text="(2)Compra\nPunt", on_release=coin_to_point, size_hint_y=None, height=BtnHeight)) args["difficulty"] = 1 args["difficultyCost"] = 5 + def raiseDifficulty(btn): - if args["coins"] >= args["difficultyCost"]*args["difficulty"]: - coins_add(args, -args["difficulty"]*args["difficultyCost"]) + if args["coins"] >= args["difficultyCost"] * args["difficulty"]: + coins_add(args, -args["difficulty"] * args["difficultyCost"]) args["difficulty"] += 1 - if args["difficulty"] == 3: args["difficultyButton"].disabled = True - else: args["difficultyButton"].text = "({})Augmenta\nDificultat".format(args["difficultyCost"] * args["difficulty"]) + if args["difficulty"] == 3: + args["difficultyButton"].disabled = True + else: + args["difficultyButton"].text = "({})Augmenta\nDificultat".format( + args["difficultyCost"] * args["difficulty"]) if args["coinSchedule"] is not None: args["coinSchedule"].cancel() args["newCoin"] = True ps.animation.frame = 1 - #args["kHarm"] += 2 - ps.QSystem.setState(mathPhysics.eigenvectorHarmonic2DGenerator(0., args["difficulty"], 0., args["difficulty"], args["kHarm"])) + # args["kHarm"] += 2 + ps.QSystem.setState(mathPhysics.eigenvectorHarmonic2DGenerator(0., args["difficulty"] - 1, 0., + args["difficulty"] - 1, + args["kHarm"])) ps.QSystem.renorm() ps.animation.manualUpdate(onlyDraw=True) ps.animation.rescalePsi() ps.QSystem.t = 0. - ps.extra_param[0:4+1] = 0. - #ps.extra_param[5] = args["kHarm"] + ps.extra_param[0:4 + 1] = 0. + # ps.extra_param[5] = args["kHarm"] - args["difficultyButton"] = Button(text="({})Augmenta\nDificultat".format(args["difficultyCost"*args["difficulty"]]), - on_release=raiseDifficulty, size_hint_y=None, height=BtnHeight) + args["difficultyButton"] = Button( + text="({})Augmenta\nDificultat".format(args["difficultyCost"] * args["difficulty"]), + on_release=raiseDifficulty, size_hint_y=None, height=BtnHeight) shop.add_widget(args["difficultyButton"]) shopBig.add_widget(shop) layout.add_widget(shopBig) - args["sizeStep"] = 1 def raiseSize(btn): - if args["coins"] >= args["sizeStep"]**2: - coins_add(args, -args["sizeStep"]**2) - args["sizeStep"]+=1 + if args["coins"] >= args["sizeStep"] ** 2: + coins_add(args, -args["sizeStep"] ** 2) + args["sizeStep"] += 1 args["goalRadius"] *= sqrt(3) args["goalCircle"].set(radius=args["goalRadius"]) - args["sizeButton"].text = "({})Objectiu\nmés gran".format(args["sizeStep"]**2) + args["sizeButton"].text = "({})Objectiu\nmés gran".format(args["sizeStep"] ** 2) - args["sizeButton"] = Button(text="({})Objectiu\nmés gran".format(args["sizeStep"]**2), + args["sizeButton"] = Button(text="({})Objectiu\nmés gran".format(args["sizeStep"] ** 2), on_release=raiseSize, size_hint_y=None, height=BtnHeight) shop.add_widget(args["sizeButton"]) @@ -1910,44 +2309,42 @@ def raiseSize(btn): return layout def extra_clean_moveGame(args): - plt.close(args["energyGraph"]) - args["energyGraph"].canvas.canvas.clear() + animate.cleanFigClose(args["energyGraph"]) if args["coinSchedule"] is not None: args["coinSchedule"].cancel() """clock = Clock.schedule_once(lambda x: print("clock"), 0.5) def unscheduleAfter(dt): clock.cancel() print("?") - Clock.schedule_once(unscheduleAfter, 1)"""# We can check like this no error is thrown for cancelling an - # already finished clock schedule - - + Clock.schedule_once(unscheduleAfter, 1)""" # We can check like this no error is thrown for cancelling an + # already finished clock schedule def moveGameInfo(): infoPopup = Popup(title="INFO", auto_dismiss=True, size_hint=(0.6, 0.6)) infoPopup.add_widget(Label(text= - "Mou el pou de potencial amb\nles fletxes del teclat ← ↑ → ↓ o WASD\n" - "Transporta així l'electró fins els cercles marcats,\non vols trobar-lo quan l'observis (en omplir-se la barra blava)\n" - "Però vigila que no el perdis! Pista: Vigila l'energia\n\nClica ràpid per atrapar les monedes i comprar millores\n" - "Quants punts pots aconseguir?", font_name="Arial")) + "Mou el pou de potencial amb\nles fletxes del teclat ← ↑ → ↓ o WASD\n" + "Transporta així l'electró fins els cercles marcats,\non vols trobar-lo quan l'observis (en omplir-se la barra blava)\n" + "Però vigila que no el perdis! Pista: Vigila l'energia\n\nClica ràpid per atrapar les monedes i comprar millores\n" + "Quants punts pots aconseguir?", font_name="Arial")) infoPopup.open() - #TextPopup("Mou el pou de potencial amb\nles fletxes del teclat ← ↑ → ↓\nTransporta l'electró per tal que es pugui\nobservar a les zones marcades").open() + # TextPopup("Mou el pou de potencial amb\nles fletxes del teclat ← ↑ → ↓\nTransporta l'electró per tal que es pugui\nobservar a les zones marcades").open() definition = { - #"QSystem": QSystem, - "initState": inicial2D, "potential":potentialHarmonicWellMoving, "extra_param":extra_param, - "drawClassical": False, "drawExpected": False, "duration": None,#10. - "extra_update": extra_update_movement, #"extraCommands": [('key_press_event', moveHarmonicWellKeyboard)], - "extraUpdates": [drawCircle], "isFocusable":False, "plotWidget":extra_info_movement, "scalePsi": False, - "showNorm": False, "showEnergy":False, "showMomentum":False, "debugTime":False, "extra_keyboard_action":extra_keyboard_movegame, - "info_action": moveGameInfo, "extra_clean":extra_clean_moveGame + # "QSystem": QSystem, + "initState": inicial2D, "potential": potentialHarmonicWellMoving, "extra_param": extra_param, + "drawClassical": False, "drawExpected": False, "duration": None, # 10. + "extra_update": extra_update_movement, # "extraCommands": [('key_press_event', moveHarmonicWellKeyboard)], + "extraUpdates": [drawCircle], "isFocusable": False, "plotWidget": extra_info_movement, "scalePsi": False, + "showNorm": False, "showEnergy": False, "showMomentum": False, "debugTime": False, + "extra_keyboard_action": extra_keyboard_movegame, + "info_action": moveGameInfo, "extra_clean": extra_clean_moveGame } """"initState": mathPhysics.eigenvectorHarmonic2DGenerator(0., 2, 0., 2, 8.), "potential": potentialHarmonicWellMoving }""" self.ids.gameSelect.add_widget( - Button(text="Transporta la partícula!", on_release=partial(self.switch, definition=definition, setExtraArgs=setMoveGame))) - + Button(text="Transporta la partícula!", + on_release=partial(self.switch, definition=definition, setExtraArgs=setMoveGame))) def switch(self, *args, definition=None, setExtraArgs=None): self.manager.get_screen("playscreen").set_self(definition, setExtraArgs=setExtraArgs) @@ -1956,11 +2353,10 @@ def switch(self, *args, definition=None, setExtraArgs=None): self.manager.get_screen("playscreen").sourceScreen = "games" - - class ColoredLabel(Label): pass + class SaveGifPopup(Popup): def __init__(self, window, duration=5., fileName="resultat", animwidth=12., animheight=7., **kwargs): self.window = window # Window holds information such as QuantumSystem and Animation @@ -1973,15 +2369,19 @@ def __init__(self, window, duration=5., fileName="resultat", animwidth=12., anim def saveAnimation(self, fName, duration, type): anim = self.window.animation animationToSave = animate.QuantumAnimation( - anim.QSystem, dtSim=anim.dtSim, stepsPerFrame=anim.stepsPerFrame, width=self.animwidth, height=self.animheight, + anim.QSystem, dtSim=anim.dtSim, stepsPerFrame=anim.stepsPerFrame, width=self.animwidth, + height=self.animheight, duration=duration, dtAnim=anim.dtAnim, callbackProgress=True, showPotential=True, varyingPotential=True, showMomentum=anim.showMomentum, showEnergy=anim.showMomentum, showNorm=anim.showNorm, scalePsi=anim.scalePsi, scaleMom=anim.scaleMom, isKivy=False, - drawClassical=anim.drawClassical, drawClassicalTrace=anim.drawClassicalTrace, drawExpected=anim.drawExpected, drawExpectedTrace=anim.drawExpectedTrace) - #animationToSave.reset_plot() - try: animationToSave.saveAnimation(fName, type) - except: TextPopup("Error, format probablement no suportat").open() + drawClassical=anim.drawClassical, drawClassicalTrace=anim.drawClassicalTrace, + drawExpected=anim.drawExpected, drawExpectedTrace=anim.drawExpectedTrace) + # animationToSave.reset_plot() + try: + animationToSave.saveAnimation(fName, type) + except: + TextPopup("Error, format probablement no suportat").open() """def on_open(self): self.ids.Nx.text = str(self.window.animation.QSystem.Nx) @@ -1995,6 +2395,7 @@ def saveAnimation(self, fName, duration, type): self.ids.dtSim.text = str(self.window.animation.dtSim)""" + class ParametersPopup(Popup): def __init__(self, window, **kwargs): self.window = window # Window holds information such as QuantumSystem and Animation @@ -2004,7 +2405,7 @@ def setPotential(self): if (self.ids.potential.on_text_validate() == 0): self.window.QSystem.changePotential(self.window.potential) self.window.animation.manualUpdate(onlyDraw=True) - #self.window.animation.reset_plot() + # self.window.animation.reset_plot() def previewPotential(self): copy = self.window.potential @@ -2024,7 +2425,7 @@ def setInitState(self): if (self.ids.initState.on_text_validate() == 0): self.window.QSystem.setState(self.window.initState) self.window.animation.manualUpdate(onlyDraw=True) - #self.window.animation.reset_plot() + # self.window.animation.reset_plot() def previewInitState(self): copy = self.window.initState @@ -2052,13 +2453,15 @@ def previewInitState(self): self.ids.dtSim.text = str(self.window.animation.dtSim)""" + class DebugPopup(Popup): def __init__(self, window, **kwargs): self.window = window # Window holds information such as QuantumSystem and Animation super(DebugPopup, self).__init__(**kwargs) + class PlayScreen(Screen): - def __init__(self, sourceScreen = "examples", **kwargs): + def __init__(self, sourceScreen="examples", **kwargs): super(PlayScreen, self).__init__(**kwargs) self.sourceScreen = sourceScreen self.extra_param = np.array([1.]) @@ -2068,14 +2471,13 @@ def __init__(self, sourceScreen = "examples", **kwargs): # Keyboard: https://kivy.org/doc/stable/api-kivy.core.window.html self.extra_keyboard_action = None self._keyboard = None - #self._keyboard = Window.request_keyboard(self._keyboard_closed, self) + # self._keyboard = Window.request_keyboard(self._keyboard_closed, self) # self._keyboard.bind(on_key_down=self._on_keyboard_down) def _keyboard_closed(self): self._keyboard.unbind(on_key_down=self._on_keyboard_down) self._keyboard = None - def _on_keyboard_down(self, keyboard, keycode, text, modifiers): if self.extra_keyboard_action is not None: self.extra_keyboard_action(self, keyboard, keycode, text, modifiers) @@ -2087,8 +2489,8 @@ def on_enter(self, *args): if self.extra_on_enter is not None: self.extra_on_enter(self.extraArgs) def set_self(self, definition=None, setExtraArgs=None): - #self._keyboard = Window.request_keyboard(self._keyboard_closed, self) - #self._keyboard.bind(on_key_down=self._on_keyboard_down) + # self._keyboard = Window.request_keyboard(self._keyboard_closed, self) + # self._keyboard.bind(on_key_down=self._on_keyboard_down) if definition == None: definition = {} @@ -2099,12 +2501,12 @@ def set_self(self, definition=None, setExtraArgs=None): self.extra_on_enter = None self.paused = True - self.Nx = 200; self.Ny = 200 + self.Nx = 200; + self.Ny = 200 L = 10. self.x0, self.y0 = -L, -L self.xf, self.yf = L, L - # We allow 16 global variables, which can be named self.nVar = 16 self.extra_param = np.zeros(self.nVar, dtype=np.float64) @@ -2132,7 +2534,7 @@ def set_self(self, definition=None, setExtraArgs=None): self.customOperator = None self.renormStep = False - + self.step = 'fastest' for key in definition: vars(self)[key] = definition[key] @@ -2143,7 +2545,8 @@ def set_self(self, definition=None, setExtraArgs=None): self.QSystem = mathPhysics.QuantumSystem2D(self.Nx, self.Ny, self.x0, self.y0, self.xf, self.yf, self.initState, potential=self.potential, extra_param=self.extra_param, - renormStep=self.renormStep, customOperator=self.customOperator) + renormStep=self.renormStep, customOperator=self.customOperator, + step=self.step) self.savedStates = [] self.tempState = {"psi": self.QSystem.psi, "x0": self.QSystem.x0, "xf": self.QSystem.xf @@ -2155,13 +2558,16 @@ def set_self(self, definition=None, setExtraArgs=None): if setExtraArgs is not None: setExtraArgs(self.extraArgs) - animKeys = {"dtSim":0.01, "dtAnim":0.04, "stepsPerFrame":0, "debugTime":False, "duration":None, - "showPotential":True, "varyingPotential":True, "showMomentum":True, "showEnergy":True, "showNorm":False, - "scalePsi":True, "scaleMom":True, "zoomMom":1.,"scalePot":True, "isFocusable":True, - "drawClassical":True, "drawClassicalTrace":False, "drawExpected":True, "drawExpectedTrace":False, - "extraCommands":[], "extraUpdates":[], "extraUpdatesStart":False, "extraUpdatesUpdate":False, - "unit_dist":unit_dist,"unit_mom":unit_mom,"unit_time":unit_time,"unit_energy":unit_energy, - "customPlot":None, "customPlotUpdate":False, "customPlotFull":False} + animKeys = {"dtSim": 0.01, "dtAnim": 0.04, "stepsPerFrame": 0, "debugTime": False, "duration": None, + "showPotential": True, "varyingPotential": True, "showMomentum": True, "showEnergy": True, + "showNorm": False, + "scalePsi": True, "scaleMom": True, "zoomMom": 1., "scalePot": True, "isFocusable": True, + "drawClassical": True, "drawClassicalTrace": False, "drawExpected": True, + "drawExpectedTrace": False, + "extraCommands": [], "extraUpdates": [], "extraUpdatesStart": False, "extraUpdatesUpdate": False, + "unit_dist": unit_dist, "unit_mom": unit_mom, "unit_time": unit_time, "unit_energy": unit_energy, + "customPlot": None, "customPlotUpdate": False, "customPlotFull": False, + "psiRepresentation": "mod2"} for key in animKeys: # Changes value to definition (dictionary), but if it's not in the dictionary leaves it as is (default) animKeys[key] = definition.get(key, animKeys[key]) @@ -2177,12 +2583,13 @@ def set_self(self, definition=None, setExtraArgs=None): # ADD HERE EXTRA THINGS LIKE SLIDERS? CAN PASS DOWN CUSTOM WIDGET ###### if "plotWidget" in definition: - if callable(definition["plotWidget"]): self.plotBox.add_widget(definition["plotWidget"](self.extraArgs, self)) # Careful here with callable() - else: self.plotBox.add_widget(definition["plotWidget"]) + if callable(definition["plotWidget"]): + self.plotBox.add_widget(definition["plotWidget"](self.extraArgs, self)) # Careful here with callable() + else: + self.plotBox.add_widget(definition["plotWidget"]) ###### - buttonBox = BoxLayout(size_hint=(1, 0.2), orientation="horizontal", padding=10, spacing=20) if self.info_action is not None: @@ -2191,11 +2598,11 @@ def set_self(self, definition=None, setExtraArgs=None): buttonBox.add_widget(resetButton) self.playButton = PlayButton(text="", state='normal' if self.paused else 'down', - on_press= lambda x: self.startPlaying() if self.paused else self.stopPlaying()) + on_press=lambda x: self.startPlaying() if self.paused else self.stopPlaying()) buttonBox.add_widget(self.playButton) - returnButton = ReturnButton(on_press = self.goBack)#, text="Retorna enrere") + returnButton = ReturnButton(on_press=self.goBack) # , text="Retorna enrere") buttonBox.add_widget(returnButton) mainBox = BoxLayout(orientation="vertical") @@ -2205,17 +2612,18 @@ def set_self(self, definition=None, setExtraArgs=None): self.add_widget(mainBox) def goBack(self, *args): - #self.clean() Clean done on_leave + # self.clean() Clean done on_leave self.manager.transition.direction = "right" self.manager.current = self.sourceScreen def clean(self): self.stopPlaying() - self.animation.fig.canvas.canvas.clear() + # self.animation.fig.canvas.canvas.clear() + animate.cleanFigClose(self.animation.fig) + self.plotRelative.clear_widgets() self.plotBox.clear_widgets() self.clear_widgets() - plt.close(self.animation.fig) if self.extra_clean is not None: self.extra_clean(self.extraArgs) ### Trying to fix memory leak, to no avail @@ -2233,11 +2641,13 @@ def resetAll(self, *args): self.clean() self.set_self(self.definition, setExtraArgs=self.setExtraArgs) Clock.schedule_once(lambda *args: Clock.schedule_once(self.on_enter)) # Everything is drawn next frame - # so we need to wait 2 frames to resize things + # so we need to wait 2 frames to resize things def stopPlaying(self): - try: self.schedule.cancel() - except: pass + try: + self.schedule.cancel() + except: + pass self.paused = True self.animation.paused = True if self._keyboard is not None: self._keyboard.unbind(on_key_down=self._on_keyboard_down) @@ -2266,12 +2676,10 @@ def setVarDict(self): {self.paramNames[i]: "extra_param[{}]".format(i) for i in range(self.nVar) if self.paramNames[i] != ""}""" - - class quantumMovementApp(App): def build(self): return WindowManager() - #return kv + # return kv if __name__ == "__main__": diff --git a/QuantumMovement/quantummovement.kv b/QuantumMovement/quantummovement.kv index c2789dba..883f3f73 100644 --- a/QuantumMovement/quantummovement.kv +++ b/QuantumMovement/quantummovement.kv @@ -16,6 +16,21 @@ # Ellipse: # pos: self.pos # size: self.size +: + canvas.before: + Color: + rgba: 240/255,240/255,240/255,255/255 + Rectangle: + pos: self.pos + size: self.size +: + Image: + fit_mode: "scale-down" + source: self.parent.image_src + center_x: self.parent.center_x + center_y: self.parent.center_y + width: self.parent.width*self.parent.size_hint_image + #allow_stretch: True : size_hint_min_y: dp(40) @@ -120,9 +135,6 @@ title: "" auto_dismiss: True size_hint: 0.5, 0.5 - BoxLayout: - id: plot - padding:20 : auto_dismiss: True @@ -140,9 +152,11 @@ id: stateName multiline: False text: "nouEstat" + on_text_validate: StateExpression: id: stateCreator holder: root + varDict: root.window.variablesDict Button: size_hint_x: 0.15 text: "Crea\nEstat" @@ -196,6 +210,7 @@ attribute: "maxiter" condition: "positive" Button: + id: eigenButton size_hint_x: 0.6 text: "Busca següent\nEigenstat" on_release: @@ -560,14 +575,26 @@ active: root.window.animation.showNorm on_active: root.window.animation.reset_plot(showNorm=self.active) - Button: - id: eigen - text: "Eigenstates" - on_release: - root.window.stopPlaying() - Factory.SavedEigenstatesPopup(root.window).open() - #root.QSystem.approximateEigenstate() - #root.animation.manualUpdate(onlyDraw=True) + GridLayout: + cols: 3 + Button: + id: eigen + text: "Eigenstates" + on_release: + root.window.stopPlaying() + Factory.SavedEigenstatesPopup(root.window).open() + #root.QSystem.approximateEigenstate() + #root.animation.manualUpdate(onlyDraw=True) + ToggleButton: + text: "Rescala psi\ncada pas" + state: 'down' if root.window.animation.scalePsi else 'normal' + on_press: + root.window.animation.reset_plot(scalePsi=not root.window.animation.scalePsi, forceRecreate=True) + ToggleButton: + text: "Rescala moment\ncada pas" + state: 'down' if root.window.animation.scaleMom else 'normal' + on_press: + root.window.animation.reset_plot(scaleMom = not root.window.animation.scaleMom, forceRecreate=True) :