VAMPyR introduction

1 Prologue

Import the packages like so

from vampyr import vampyr3d as vp
import matplotlib.pyplot as plt
import numpy as np

If you want to find out how to use a function, or what we have in the package use the function below

help(vp.ScalingProjector)

In case the help was not enough, ask us for help. For this you need only your voice.

2 Math in VAMPyR

The following exercises are introductions in working with functions, application of operators and some equation solving.

2.1 Projecting a function

In this task you will create a function, project it at different precisions and basis orders, and plot it against its projection. This is meant to give you understanding in how one can approach completeness in multiwavelets.

First we start by initializing the multiresolution analysis (mra) by using the following command

mra = vp.MultiResolutionAnalysis(box=[-20,20], order=7)

Here you will change the order value later. This defines the polynomial order of the basis.

An allowed function for projection in vampyr needs to have a single coordinate argument and it has to output a float value. You will construct a gaussian function defined by $\( f(\mathbf{r}) = \beta \exp(-\alpha[r_x^2 + r_y^2 + r_z^2]) \)$ the alpha value is a measure of how sharp the gaussian function is, and the beta is just a normalization factor, for this exercise you can set it to anything. Fill out the following python function for this (hint: you can use np.exp).

def f(r):
    # implement your function here (x,y,z is accessed as r[0],r[1],r[2])
    return

You should then be able to plot this function by running the next cell

r_x = np.linspace(-0.99, 0.99, 1000) # create an evenly spaced set of points between -0.99 and 0.99
r_y = np.zeros(1000)
r_z = np.zeros(1000)
r = [r_x, r_y, r_z]

plt.plot(r_x, f(r))

Now we want to project this function, the way we do that is by creating a projection operator vp.ScalingProjector and apply it on the function. We apply this operator by using the function as a function argument of it, P(f) for example.

Try creating one and project the function you created earlier.

hint: use the help function to fill out the blanks)

### write your code here
P_eps = 

f_tree =    # <-- here is where you will put the projected function

when you are done you can run the following cell to plot your projected function.

f_plt = [ f_tree([x, 0.0, 0.0]) for x in r_x ]
plt.plot(r_x, f(r), "tab:blue") # your function
plt.plot(r_x, f_plt, "tab:red") # your projected function
plt.show()

Now that you have projected the function, try changing three values, the alpha value of your original function, the order of the basis and the precision (prec) at which you project the function. What you should see is that the sharper the function, the higher order and precision you need to represent it properly.

One last thing I need to point out before we go to the next exercise is that it is convention to call our projected functions trees, as they all are vp.FunctionTree type.

2.2 Arithmethic

Arithmethic in vampyr is mostly the same as in normal python. Lets work with the previous function you created.

First you can normalize your projected function by calling the function normalize() of f_tree. Try it in the next cell

# normalize your function here

Before we go on into the rest of the exercise, it is important to note that whenever a function is shown, and it doesn’t start with vp., it is assumed that it must be called as a member function of the tree, that is f_tree.normalize(), otherwise all functions must be called by using the vp. prefix.

Start by multiplying it by a factor of \(2\) and assigning it to a new tree g_tree

g_tree =  # <--- assign your rescaled function to this

You should be able to compare it against f_tree by runnning the following cell

g_plt = [ g_tree([x, 0.0, 0.0]) for x in r_x ]
plt.plot(r_x, g_plt, "tab:red")
plt.plot(r_x, f_plt, "tab:blue")
plt.show()

try adding f_tree and g_tree together

h_tree =  # <--- assign your added functions to this

again, you can test this by plotting it against the two other functions by running the following cell

h_plt = [ h_tree([x, 0.0, 0.0]) for x in r_x ]

plt.plot(r_x, h_plt, "tab:orange")
plt.plot(r_x, g_plt, "tab:red")
plt.plot(r_x, f_plt, "tab:blue")
plt.show()

You can also check that the norm of the different trees makes sense. This is done by calling the norm() value of your tree, try printing it below using the print function. The f_tree should be normalized, so its norm should be about \(1.0\)

print() # print the norm of your functions

Try multiplying your trees with the * operator, it should return a new tree which is the product of both, and you can check its norm again.

i_tree =  # <----- assign you multiplied functions to this

Finally, you can take the dot product, also called inner product, between two trees by using the vp.dot function. You can try this in the cell below. print out the norm of the product as shown above, is it what you would expect?

# take the dot product of f_tree and g_tree

2.3 Operators

You have already worked with the vp.ScalingProjector and will work with the vp.HelmholtzOperator. In this part we will work with the Derivative operator, defined as vp.ABGVDerivative (named after the authors of the paper: Alpert, Beylkin, Gines and Vozovoi). Try constructing (with the help of help()) and applying the operator to your f_tree

D =  # <--- initialize the ABGVDerivative operator here
df_tree =  #  <-- set your differentiated f_tree here

Run the cell below to plot your differentiated function

df_plt = [ df_tree([x, 0.0, 0.0]) for x in r_x ]
plt.plot(r_x, df_plt, "tab:blue") # your differentiated function
plt.plot(r_x, f_plt, "tab:red") # your function tree
plt.show()