wavefunction_cartesian_effective_mass.c++

Go to the documentation of this file.

/*
 *Provides a nice class to house an effective mass wavefunction.

    Copyright (C) 2010 Joseph Pingenot

    This program is free software: you can redistribute it and/or modify
    it under the terms of the GNU Affero General Public License as published by
    the Free Software Foundation, either version 3 of the License, or
    (at your option) any later version.

    This program is distributed in the hope that it will be useful,
    but WITHOUT ANY WARRANTY; without even the implied warranty of
    MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
    GNU Affero General Public License for more details.

    You should have received a copy of the GNU Affero General Public License
    along with this program.  If not, see <http://www.gnu.org/licenses/>.

*/

#include <libkdotp/wavefunction_cartesian_effective_mass.h++>
#include <math.h>

namespace kdotp {
  namespace libkdotp {
    namespace wavefunction {

      #define WAVEFUNCTION_CARTESIAN_EFFECTIVE_MASS_GERROR_DOMAIN 29
      

      cartesian_effective_mass::cartesian_effective_mass(double* eigenvector_block, unsigned int eigvec, double dx, unsigned int L, double energy) : E(energy), dx(dx), L(L) {
  wf = new double[L];
  for(unsigned int i=0; i<L; i++) wf[i] = eigenvector_block[L*eigvec+i];
      }

      /*Although this is authoritative, it's probably faster to use L*eigvec+i directly, so the L*eigvec is pulled out of the loop*/
      ___const unsigned int get_cartesian_effective_mass_block_index(unsigned int i, unsigned int eigvec, unsigned int L) {
  return L*eigvec + i;
      }

      double get_average_with_cartesian_effective_mass_block(double* block, unsigned int eigvec, unsigned int L, double* func) {
  double avg = 0;
  for(unsigned int i=0; i<L; i++) avg += block[L*eigvec + i]*block[L*eigvec + i]*func[i];
  return avg;
      }

      /**Sets the points to the density from this wavefunction.
       *\param block the block of eigenvalues
       *\param eigvec which eigenvector to use in the block
       *\param L the number of points in each wave function in the block
       *\param dx inter-node spacing
       *\param density pre-allocated structure to store density in.
       *\return pointer to GError 
       *\note the density has already been allocated!
       */
      GError* init_density_1d_from_cartesian_effective_mass_block(double* block, unsigned int eigvec, unsigned int L, double dx, struct density_1d *density) {
  if(unlikely(density->N != L)) return g_error_new(WAVEFUNCTION_CARTESIAN_EFFECTIVE_MASS_GERROR_DOMAIN, 1, "::kdotp::libkdotp::wavefunction::wavefunction_cartesian_effective_mass.c++::init_density_1d_from_carteisan_effective_mass_block: allocated size of the density's storage (%u) and the number of points in the %u-th eigenvector (%u) don't agree!", density->N, eigvec, L);
  density->dx = dx;
  density->N = L;
  for(unsigned int i=0; i<L; i++) density->storage[i] = block[L*eigvec + i]*block[L*eigvec + i];
  return NULL;
      }

      double get_cartesian_effective_mass_block_max_delta(double* block1, unsigned vec1, double* block2, unsigned vec2, unsigned int L) {
  //fprintf(stderr, "In get_cartesian_effective_mass_block_max_delta:\n");
  if((block1 == NULL) || (block2 == NULL)) return -1;
  double max_change = fabs(block1[L*vec1]*block1[L*vec1] - block2[L*vec2]*block2[L*vec2]);
  register double change;
  for(unsigned int i=1; i<L; i++) {
    change = fabs(block1[L*vec1+i]*block1[L*vec1+i] - block2[L*vec2+i]*block2[L*vec2+i]);
    //fprintf(stderr, "\t%u: change=%g maxchange=%g (%g-%g)\n", i, change, max_change, block1[L*vec1+i]*block1[L*vec1+i], block2[L*vec2+i]*block2[L*vec2+i]);
    if(change > max_change) max_change = change;
  }
  //fprintf(stderr, "\tmaxchange=%g\n",  max_change);
  return fabs(max_change);
      }
    }
  }
}