Least Action

Nontrivializing triviality..and vice versa.

Archive for the ‘Physics’ Category

Fun with Homology

leave a comment »

There is a theorem (currently being attributed to Wikipedia, but I’m sure I can do better given more time) which states that

All closed surfaces can be produced by gluing the sides of some polygon and all even-sided polygons (2n-gons) can be glued to make different manifolds.

Conversely, a closed surface with n non-zero classes can be cut into a 2n-gon.

Two interesting cases of this are:

  1. Gluing opposite sides of a hexagon produces a torus T^2.
  2. Gluing opposite sides of an octagon produces a surface with two holes, topologically equivalent to a torus with two holes.

I had trouble visualizing this on a piece of paper, so I found two videos which are fascinating and instructive, respectively.

The two-torus from a hexagon

The genus-2 Riemann surface from an octagon

I would like to figure out how one can make such animations, and generalizations of these, using Mathematica or Sagemath.

There are a bunch of other very cool examples on the Youtube channels of these users. Kudos to them for making such instructive videos!

PS – I see that $\LaTeX$ on WordPress has become (or is still?) very sloppy! ūüė¶


Written by Vivek

October 20, 2016 at 21:57

Divergence Theorem in Complex Coordinates

leave a comment »

The divergence theorem in complex coordinates,

\int_R d^2{z} (\partial_z v^z + \partial_{\bar{z}}v^{\bar{z}}) = i \oint_{\partial R}(v^z d\bar{z} - v^{\bar{z}}dz)

(where the contour integral circles the region R counterclockwise) appears in the context of two dimensional conformal field theory, to derive Noether’s Theorem and the Ward Identity for a conformally invariant scalar field theory (for example), and is useful in general in 2D CFT/string theory. This is equation (2.1.9) of Polchinski’s volume 1, but a proof is not given in the book.

This is straightforward to prove by converting both sides separately to Cartesian coordinates (\sigma^1, \sigma^2), through

z = \sigma^1 + i \sigma^2
\bar{z} = \sigma^1 - i \sigma^2

\partial_z = \frac{1}{2}(\partial_1 - i \partial_2)

\partial_{\bar{z}} = \frac{1}{2}(\partial_1 + i \partial_2)

d^2 z = 2 d\sigma^1 d\sigma^2 = 2 d^2 \sigma

and using the Green’s theorem in the plane

\oint_{\partial R}(L dx + M dy) = \int \int_{R} \left(\frac{\partial M}{\partial x} - \frac{\partial L}{\partial y}\right) dx dy

with the identifications

x \rightarrow \sigma^1, y \rightarrow \sigma^2
L \rightarrow -v^2, M \rightarrow v^1

There is perhaps a faster and more elegant way of doing this directly in the complex plane, but this particular line of reasoning makes contact with the underlying Green’s theorem in the plane, which is more familiar from real analysis.

The Yang Mills Field Strength: a motivation for differential forms

leave a comment »

Although their importance and applications are widely recognized in theoretical physics, differential forms are not part of the standard curriculum in physics courses, except for the rare mention in a general relativity course, e.g. if you have used the book by Sean Carroll. In this blog post, we describe the construction of the field strength tensor for Yang Mills Theory, using differential forms. The derivation (which is based on the one given by Tony Zee in his QFT book) is elegant, and mathematically less tedious than the more conventional derivation based on the matrix gauge field. It also illustrates the power of the differential forms approach.

To begin with, let me define a normalized matrix gauge field A_\mu which is -i times the usual matrix gauge field. So in the equations that appear below, the covariant derivative has no lurking factors of -i. This greatly simplifies the algebra for us, as we don’t have to keep track of conjugations and sign changes associated with i. The covariant derivative is thus D_\mu = \partial_\mu + A_\mu in terms of this new gauge field. The matrix 1-form A is

A = A_\mu dx^\mu

So, A^2 = A_\mu A_\nu dx^\mu dx^\nu. But since dx^\mu dx^\nu = -dx^\nu dx^\mu, only the antisymmetric part of the product survives and hence we can write

A^2 = \frac{1}{2}[A_\mu, A_\nu] dx^\mu dx^\nu

We want to construct an appropriate 2-form F = \frac{1}{2}F_{\mu\nu}dx^\mu dx^\nu from this 1-form A. Now, if d denotes the exterior derivative, then dA is a 2-form, as is A^2. These are the only two forms we can construct from A. So, F must be a linear combination of these two forms. This is a very simple, and neat argument!

Now, the transformation law for the gauge potential is

A \rightarrow U A U^\dagger + U dU^\dagger

where U is a 0 form (so dU^\dagger = \partial_\mu U^\dagger dx^\mu). Applying d to the transformation law, we get

dA \rightarrow U dA U^\dagger + dU A U^\dagger - U A dU^\dagger + dU dU^\dagger

where the negative sign in the third term comes from moving the 1-form d past the 1-form A. Squaring the transformation law yields

A^2 \rightarrow UA^2 U^\dagger + U A dU ^\dagger + U dU^\dagger U A U^\dagger + U dU^\dagger U dU^\dagger

Now, $UU^\dagger = 1$, so applying d again to both sides we get U dU^\dagger = -dU U^\dagger. So, we can write the square transformation law as

A^2 \rightarrow U A^2 U^\dagger + U A dU^\dagger - dU A U^\dagger - dU dU^\dagger
whereas if we recall the expression for the transformation of dA, it was just
dA \rightarrow U dA U^\dagger + dU A U^\dagger - U A dU^\dagger + dU dU^\dagger

Clearly if we merely add A^2 and dA, the last 3 terms on the RHS of each cancel out, and we get

A^2 + dA \rightarrow U(A^2 + dA)U^\dagger

which is the expected transformation law for a field strength of the form F = A^2 + dA:

F \rightarrow U F U^{\dagger}

The differential form approach uses compact notation that suppresses the Lorentz index \mu as well as the group index a, and gives us a fleeting glimpse into the connection between gauge theory and fibre bundles.

For a gentle yet semi-rigorous introduction to differential forms, the reader is referred to the book on General Relativity by Sean Carroll.

Written by Vivek

June 11, 2014 at 10:53

Gnuplot with C/C++

leave a comment »

Sometimes, it is convenient to call gnuplot from within your own C/C++ code rather than having the code create a data file, and manually execute gnuplot every time to re-plot using the contents of the file. This is easy if you use Linux: just use a pipe. For instance:

#include <iostream>
#include <stdio.h>
#include <cstdlib>

using namespace std;

int main(){
// Code for gnuplot pipe
FILE *pipe = popen("gnuplot -persist", "w");
fprintf(pipe, "\n");

// Your code goes here

// Don't forget to close the pipe
return 0;

Note that we have had to use C routines to call gnuplot. But there’s a problem with this code: gnuplot will generate the plot only when the pipe is closed. So if you have multiple plot statements, or if you want to refresh the same plot (as for instance, in an animation scenario) you will find that the above approach will result in the plots getting asynchronously updated right at the end.

So how does one ensure that after every plot command, a plot is actually generated (or an existing one refreshed) Immediately? The key is to flush the buffer using fflush. Simply insert


below every plot call to gnuplot. If you have an animation, then in order to make the transitions smoother, you might want to include a delay subroutine. A simple one can be built using the clock routines in time.h. For an example, click here.

Thanks to Vidushi Sharma for bringing these issues to my attention.

Written by Vivek

June 22, 2012 at 14:13

Why V(x) = -1/x^2 has no bound state

leave a comment »

Here we present a scaling-based argument to show that the attractive potential

V(x) = -\frac{\lambda}{x^2}

(\lambda > 0), has no bound states (i.e. states with energy E < 0). Consider the Time Independent Schrodinger equation for this potential, which is the eigenvalue equation for the corresponding Hamiltonian,

-\frac{\hbar^2}{2m}\frac{d^2\psi}{dx^2} - \frac{\lambda}{x^2}V(x)\psi = E\psi

This can be rearranged as

\frac{d^2\psi}{dx^2} + \frac{2m\lambda}{\hbar^2}\frac{1}{x^2}\psi = -\frac{2mE}{\hbar^2}\psi

Now, it is easy to see that the quantity


is dimensionless. So, this problem has no independent scale, even though naively one might think that \lambda specifies a scale for this problem. We claim that for such a system, there can be no bound state. This is proved below.

Suppose we perform the scale transformation x \rightarrow \alpha x where \alpha is some nonzero real number, we see from the equation above that if E is an eigenvalue, then so is \alpha^2 E.

Suppose now that a bound ground state exists, with energy E_{G}. By definition E_{G} < 0. Then scale invariance implies that \alpha^2 E_{G} must also be the energy of some bound state. But

\alpha^2 E_{G} < E_{G}

as multiplying the negative ground state energy by a positive number only makes it more negative. This contradicts the fact that the ground state has energy E_{G}. In fact, we can make a stronger statement, viz. the ground state energy can be made as small as we want. Therefore, there is no finite energy ground state for this system, and consequently there can be no bound states.

Note that there is no such problem with scattering states, i.e. states with positive energy. One can take an arbitrarily small positive energy scattering state and from it obtain valid energies of the continuum of higher energy scattering states by performing a scale transformation.

Incidentally, this is why potentials like -\lambda[\delta(x)]^2 and -\lambda\delta(x)/x also have no bound states.

Written by Vivek

February 8, 2012 at 08:48

Evaluation of Complex Error Functions erf(z) using GSL/C

with 5 comments

For a small QFT calculation, I needed to numerically evaluate the imaginary error function erfi(x) = erf(i x). But it turns out that GSL (and most other numerical recipe code I could find) can only deal with erf(x), where x is real. Here’s a poor man’s implementation of erf(z) through a standard Taylor expansion,

erf(z) = \frac{2}{\sqrt{\pi}}\sum_{n = 0}^{\infty}\frac{(-1)^n z^{2n+1}}{n!(2n+1)}

The catch here is to deal with propagation of errors in the complex Taylor series, and to also somehow benchmark the results. Well, I haven’t been able to think about this yet, but I was able to confirm that for real arguments, my Taylor series code is about as good as the GSL error function gsl_sf_erf(x).

I am working on a CUDA implementation of this now, because in my project, I need to perform a numerical integration over the error function, which is quite intensive even for Mathematica. For now, I’m just sharing the serial implementation. (Hint: if you can write wrapper for each of the gsl functions inside my Taylor series calculating method, you can embed them within a __global__ kernel call through a struct. Alternatively — and less fun — you can just parallelize the for loop.)

/* Function to compute erf(z) using a Taylor series expansion
/* Author: Vivek Saxena
/* Last updated: January 28, 2011 21:11 hrs */
#include <gsl/gsl_complex_math.h>
#include <gsl/gsl_sf_erf.h>
#include <gsl/gsl_sf_pow_int.h>
#include <gsl/gsl_sf_gamma.h>
#include <stdio.h>

#define PI 3.1415926543
double cz = 2/sqrt(PI);
const int TERMS = 10; // no of terms to use in the Taylor series

gsl_complex erfTaylor( gsl_complex z, int trunc ){
  gsl_complex res = gsl_complex_rect(0,0),
  num  = gsl_complex_rect(0,0),
  den  = gsl_complex_rect(1,0),
  snum = gsl_complex_rect(1,0),
  temp = gsl_complex_rect(0,0);
  for(int i = 0; i < trunc; i++){
     snum = gsl_complex_rect( cz * gsl_sf_pow_int(-1, i), 0 );
     num  = gsl_complex_mul(snum, gsl_complex_pow_real(z, 2*i+1));
     den  = gsl_complex_rect((2*i + 1)*gsl_sf_fact(i),0);
     temp = gsl_complex_div(num, den);
     res  = gsl_complex_add(res, temp);

 return res;

int main ( void ){
   printf( "Real error function\n\n");
   for ( float i = 0; i <= 1; i += 0.01 ){
     float gslerror = gsl_sf_erf(i);
     float taylor   = GSL_REAL( erfTaylor( gsl_complex_rect(i, 0), 10 ) );
     printf("erf(%f): gsl = %f, taylor = %f, mag error = %f\n", i, gslerror, taylor, abs(gslerror-taylor));

   gsl_complex t, arg;
   printf( "\n\nImaginary error function\n\n");
   for (float i = 0; i <= 1; i += 0.01 ){
      /* this would be your generic argument z in erf(z).
       * so if z = x + iy, then
       * arg = gsl_complex_rect(x, y);
      arg = gsl_complex_rect(0, i);
      t   = erfTaylor( arg, TERMS );
      printf("erf(%f + i %f) = %f + i %f\n", GSL_REAL(arg), GSL_IMAG(arg), GSL_REAL(t), GSL_IMAG(t));
   return 0;

Written by Vivek

January 28, 2011 at 21:12

A Mathematica notebook for Angular Momentum Matrix Algebra

leave a comment »

This is a Mathematica notebook for angular momentum matrix algebra. I have not yet added a component to add angular momenta, which I plan to do soon.

Feel free to use and distribute the notebook. Please just add a citation and a link to this webpage (also optional, of course).

Angular Momentum Matrix Algebra Notebook

Written by Vivek

November 17, 2010 at 19:12