V1: VAMPyR introduction#
Prologue#
The VAMPyR package can be imported in the following way, along with any other useful Python package you want to use. Remember to choose the appropriate dimensionality in the vampyr package name.
from vampyr import vampyr3d as vp
import matplotlib.pyplot as plt
import numpy as np
Some documentation is available through the help
function, e.g.
help(vp.ScalingProjector)
MultiResolution Analysis (MRA)#
The MRA defines the computational domain and the polynomial basis. The domain needs to be either symmetric around the origin \([-L,L]^D\) or in the positive hyperoctant \([0,L]^D\).
from vampyr import vampyr3d as vp
import matplotlib.pyplot as plt
import numpy as np
k = 5 # Polynomial order 0 < k < 40
L = 1 # Box size needs to be [-L,L] or [0,L]
epsilon = 1.0e-3 # Precision parameter
# Define MRA in the box [-L,L]^3
MRA1 = vp.MultiResolutionAnalysis(box=[-L,L], order=k)
print(MRA1)
# Define MRA in the box [0,L]^3
MRA2 = vp.MultiResolutionAnalysis(box=[0,L], order=k)
print(MRA2)
Function projection#
In order to represent functions in the multiwavelet basis we define a ScalingProjector
on the MRA. The projector can either be given a fixed uniform scale
or an adaptive prec
parameter.
Example: Fixed scale projection 1D#
from vampyr import vampyr1d as vp
import matplotlib.pyplot as plt
import numpy as np
k = 3 # Polynomial order
L = 1 # Box size [-L,L]
n = 3 # Scale parameter for uniform projection 2^-n
# Define fixed scale projector onto MRA
MRA = vp.MultiResolutionAnalysis(box=[-L,L], order=k)
P_n = vp.ScalingProjector(MRA, scale=n)
# Define analytic function
def f(r):
x = r[0] # r is a coordinate vector, even in 1D
return np.exp(-5*x**2)*np.cos(8*np.pi*x)
# Project function onto MRA
f_n = P_n(f)
# Plot function in the interval x=[-L,L]
r_x = np.linspace(-L, L, 1000)
plt_ana = [f([x]) for x in r_x]
plt_mra = [f_n([x]) for x in r_x]
plt.plot(r_x, plt_ana, "tab:blue") # analytic function
plt.plot(r_x, plt_mra, "tab:red") # projected function
plt.show()
Example: Adaptive projection 3D#
from vampyr import vampyr3d as vp
import matplotlib.pyplot as plt
import numpy as np
k = 3 # Polynomial order
L = 2 # Box size [0,L]
epsilon = 1.0e-3 # Precision parameter
# Define adaptive projector onto MRA
MRA = vp.MultiResolutionAnalysis(box=[0,L], order=k)
P_eps = vp.ScalingProjector(MRA, prec=epsilon)
# Define analytic function
def f(r):
R = np.sqrt(r[0]**2 + r[1]**2 + r[2]**2)
return np.exp(-R**2)*np.sin(4*np.pi*R)
# Project function onto MRA
f_eps = P_eps(f)
# Plot function in the interval x=[0,L]
r_x = np.linspace(0, L, 1000)
plt_ana = [f([x, 0.0, 0.0]) for x in r_x]
plt_mra = [f_eps([x, 0.0, 0.0]) for x in r_x]
plt.plot(r_x, plt_ana, "tab:blue") # analytic function
plt.plot(r_x, plt_mra, "tab:red") # projected function
plt.show()
Arithmethic operations#
Arithmethic in vampyr is mostly the same as in normal python. First we define two functions to be used in the following examples.
from vampyr import vampyr1d as vp
import matplotlib.pyplot as plt
import numpy as np
k = 5 # Polynomial order
L = 4 # Box size [-L,L]
epsilon = 1.0e-3 # Precision parameter
# Prepare simple plotting function for visualizations
def line_plot_1d(f, color):
r_x = np.linspace(-L, L, 1000)
plt_f = [f([x]) for x in r_x]
plt.plot(r_x, plt_f, color)
return plt
# Define fixed scale projector onto MRA
MRA = vp.MultiResolutionAnalysis(box=[-L,L], order=k)
P_eps = vp.ScalingProjector(MRA, prec=epsilon)
# Project analytic functions onto MRA
f = P_eps(lambda r : np.exp(-(r[0] - 1.0)**2))
g = P_eps(lambda r : np.exp(-(r[0] + 1.0)**2))
# Plot functions
plt = line_plot_1d(f, "tab:blue")
plt = line_plot_1d(g, "tab:red")
plt.show()
Example: scaling and normalization#
f *= np.pi # Scale f by factor pi
g.normalize() # Normalize g
# Plot functions
plt = line_plot_1d(f, "tab:blue")
plt = line_plot_1d(g, "tab:red")
plt.show()
Example: addition and multiplication#
h1 = 0.5*f - 2*g # linear combination
h2 = f * g # multiplication
h3 = g**2 # power
# Plot functions
plt = line_plot_1d(h1, "tab:blue")
plt = line_plot_1d(h2, "tab:red")
plt = line_plot_1d(h3, "tab:green")
plt.show()
Example: norms, integrals and dot products#
f_norm = f.norm()
f_sqnorm = f.squaredNorm()
f_integral = f.integrate()
f2_integral = (f**2).integrate()
f_dot_f = vp.dot(f,f)
f_dot_g = vp.dot(f,g)
print("Norm of f : {}".format(f_norm))
print("Squared norm of f : {}".format(f_sqnorm))
print("Integral of f : {}".format(f_integral))
print("Integral of f**2 : {}".format(f2_integral))
print("Dot product <f|f> : {}".format(f_dot_f))
print("Dot product <f|g> : {}".format(f_dot_g))
Operators#
Example: convolution operator#
from vampyr import vampyr3d as vp
import matplotlib.pyplot as plt
import numpy as np
k = 5 # Polynomial order
L = 10 # Box size [-L,L]
epsilon = 1.0e-3 # Precision parameter
# Prepare simple plotting function for visualizations
def line_plot_3d(f, color):
r_x = np.linspace(-L, L, 1000)
plt_f = [f([x, 0, 0]) for x in r_x]
plt.plot(r_x, plt_f, color)
return plt
# Define fixed scale projector onto MRA
MRA = vp.MultiResolutionAnalysis(box=[-L,L], order=k)
P_eps = vp.ScalingProjector(MRA, prec=epsilon)
# Project analytic functions onto MRA
rho = P_eps(lambda r : np.exp(-(r[0]**2 + r[1]**2 + r[2]**2)))
# Construct and apply Poisson operator
poisson = vp.PoissonOperator(MRA, prec=epsilon)
V = poisson(4*np.pi*rho)
# Compute electrostatic energy
E = vp.dot(V, rho)
print("Electrostatic energy : {}".format(E))
# Plot functions
plt = line_plot_3d(rho, "tab:blue")
plt = line_plot_3d(V, "tab:red")
plt.show()
Example: derivative operator#
from vampyr import vampyr1d as vp
import matplotlib.pyplot as plt
import numpy as np
k = 5 # Polynomial order
L = 4 # Box size [-L,L]
epsilon = 1.0e-3 # Precision parameter
# Prepare simple plotting function for visualizations
def line_plot_1d(f, color):
r_x = np.linspace(-L, L, 1000)
plt_f = [f([x]) for x in r_x]
plt.plot(r_x, plt_f, color)
return plt
# Define fixed scale projector onto MRA
MRA = vp.MultiResolutionAnalysis(box=[-L,L], order=k)
P_eps = vp.ScalingProjector(MRA, prec=epsilon)
# Project analytic functions onto MRA
f = P_eps(lambda r : np.exp(-r[0]**2))
# Construct and apply derivative operator. ABGV comes from the authors
# of the paper (Alpert, Beylkin, Gines, Vozovoi) and a,b are the
# boundary parameters from the same paper
D = vp.ABGVDerivative(MRA, a=0.5, b=0.5)
# Compute first and second derivative
df = D(f)
ddf = D(df)
# Plot functions
plt = line_plot_1d(f, "tab:blue")
plt = line_plot_1d(df, "tab:red")
plt = line_plot_1d(ddf, "tab:green")
plt.show()