Theory:

1. How can you find out whether GSL is installed on your box? Hint: gsl-config.

2. How can you find out which C-flags you need to use with your GSL installation?

3. How can you find out which libraries you need to link with your GSL installation?

4. Explain the syntax `command` in bash. Hint: Grave accent in programming → Unix shells; Bash command substitution.

5. Explain the syntax $(command) in bash. Hint: POSIX shell command substitution; Bash command substitution.

6. Why do you need double-dollar, $$(command), in the Makefile? Hint: GNU make: variables in recipes.

7. What is "recipe echoing"?

8. What will the following Makefile print?

   
   pwd = a string
   test:
   	@echo pwd
   	@echo `pwd`
   	@echo $(pwd)
   	@echo $$(pwd)
   

Practice: Quantum particle in a box or Standing wave in a string.

• Introduction:

  The wavefunction ψ(ξ) of a particle in a box satisfies the Schrödinger equation

  -ψ''=εψ , or Hψ=εψ, where H=-d²/dξ².

  where the prime denotes the derivative over the dimensionless variable ξ=x/L, x is the coordinate of the particle and L is the size of the box, ε is the dimensionless energy. The actual energy is E=(ℏ²ε)/(2mL²).

  The boundary condition is ψ(0)=ψ(1)=0.

  Numerically one can represent the wave-function on a grid

  ξ_i = ih, h=1/N, i=0...N,

  where the wave-function is represented by a vector

  ψ_i = ψ(x_i), i=0...N.

  The second derivative operation can then be approximated using the three point finite difference formula,

  ψ''_i = (ψ_i-1-2ψ_i+ψ_i+1)/h² .

  Now the Schrödinger equation becomes a matrix eigenvalue equation

  ∑_j=1,N-1 H_ij ψ_j = ε ψ_i, i=1...N-1,

  where—since the wave-function at the end-points is zero—the summation goes over 1...N-1 and vanishing components ψ₀ and ψ_N have been omitted.

  In matrix notation this is simply an eigenvalue equation, Hψ=εψ, where the matrix H in the three-point approximation is given as (taking into account the zero boundary condition and omitting the first and the last component of the wave-function vector)

  	(	-2	1	0	0	0	)
  	(	1	-2	1	0	0	)
  	(	0	1	-2	1	0	)
  h²H = -	(	...	...	...	...	...	)
  	(	0	0	1	-2	1	)
  	(	0	0	0	1	-2	)

  Solving the matrix eigenvalue equation amounts to matrix diagonalization.

• Exercise:
  1. Build the Hamiltonian matrix h²H as a gsl_matrix structure. Hint:

     int n=20;
     gsl_matrix *H = gsl_matrix_calloc(n,n);
     for(int i=0;i<n-1;i++){
       gsl_matrix_set(H,i,i,-2);
       gsl_matrix_set(H,i,i+1,1);
       gsl_matrix_set(H,i+1,i,1);
       }
     gsl_matrix_set(H,n-1,n-1,-2);
     gsl_matrix_scale(H,-1);
     
     

  2. Diagonalize the matrix using GSL subroutines for diagonalization of real symmetric matrices. Hint:

     gsl_eigen_symmv_workspace* w =  gsl_eigen_symmv_alloc (n);
     gsl_vector* eval = gsl_vector_alloc(n);
     gsl_matrix* evec = gsl_matrix_calloc(n,n);
     gsl_eigen_symmv(H,eval,evec,w);
     

  3. Sort the eigenvalues and eigenvectors.

     gsl_eigen_symmv_sort(eval,evec,GSL_EIGEN_SORT_VAL_ASC);
     

  4. Check that your energies are correct. Hint: the ratio should be close to unity (why?):

       for (int k=0; k < n/3; k++){ /* print first n/3 eigenvalues (rescaled) */
         double exact = pi*pi*(k+1)*(k+1)/(n+1)/(n+1);
         double ratio = gsl_vector_get(eval,k)/exact;
         fprintf (stderr, "n = %i E_n/E_exact = %g\n", k, ratio);
       }
     

  5. Plot the several lowest eigenfunction (and compare with the analytical result).

     for(int k=0;k<3;k++){
       printf("%i %f\n",0,0.0);
       for(int i=0;i<n;i++) printf("%i %f\n",i+1,gsl_matrix_get(evec,i,k));
       printf("%i %f\n",n+1,0.0);
       printf("\n\n");
       }
     

  6. You may try to investigate how the accuracy changes with the number of grid points.
  7. You may try impove the accuracy of the method by using the five-point formula for the finite difference approximation to second derivative from here.