/*
DFT++ is a density functional package developed by the research group
of Professor Tomas Arias
Copyright 1996-2003 Sohrab Ismail-Beigi
This file is part of DFT++.
DFT++ is free software; you can redistribute it and/or modify
it under the terms of the GNU General Public License as published by
the Free Software Foundation; either version 2 of the License, or
(at your option) any later version.
DFT++ 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 General Public License for more details.
You should have received a copy of the GNU General Public License
along with DFT++; if not, write to the Free Software
Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
Please see the file CREDITS for a list of authors.
For academic users, we request that publications using results obtained with
this software reference
"New algebraic formulation of density functional calculation," by Sohrab Ismail-Beigi
and T.A. Arias, Computer Physics Communications 128:1-2, 1-45 (June 2000).
and, if using the wavelet basis, further reference
"Multiresolution analysis of electronic structure: semicardinal and wavelet bases,"
T.A. Arias, Reviews of Modern Physics 71:1, 267-311 (January 1999).
and
"Robust ab initio calculation of condensed matter: transparent convergence through
semicardinal multiresolution analysis,'' I.P. Daykov, T.A. Arias, and
Torkel D. Engeness, Physical Review Letters, 90:21, 216402 (May 2003).
For your convenience, preprints of the above articles may be obtained from
http://arXiv.org/abs/cond-mat/9909130, 9805262, and 0204411, respectively.
*/
/*
* Gabor Csanyi 8/1/2001
*
*
* the rArray class holds a simple array of real numbers
*
*/
#include "header.h"
RealArray::RealArray(int length):
ndata(length)
{
}
void
RealArray::init(int length)
{
ndata = length;
if(ndata > 0){
d = (real *)malloc(ndata*sizeof(real));
//dft_log("Allocated RealArray\n");
}
else
die("Can't initiliaze RealArray.init() with length < 0\n");
}
//////////////////////////////////////
// //
// Operator functions //
// //
//////////////////////////////////////
/* Assignment: nonstandard in that it returns void. To make it standard,
* replace void -> RealArray and uncomment the return *this; */
void
RealArray::operator=(const RealArray &rd)
{
dft_log("RealArray::operator=\n");
if (ndata != rd.ndata)
die("Size mismatch in RealArray::operator=()\n");
for (int i=0; i < ndata; i++)
d[i] = rd.d[i];
/* return *this; */
}
// Assignment of scalar: all entries set to r
inline void
RealArray::operator=(real r)
{
dft_log("RealArray::operator=(real r)\n");
for (int i=0; i < ndata; i++)
d[i] = r;
}
/* Add two RealArray */
RealArray
RealArray::operator+(const RealArray &rd)
{
dft_log("RealArray::operator+\n");
if (ndata != rd.ndata)
die("Size mismatch in RealArray::operator+()\n");
RealArray rd2(*this);
for(int i = 0; i < ndata; i++)
rd2.d[i] += rd.d[i];
return rd2;
}
/* Accumulate sum of RealArray */
void
RealArray::operator+=(const RealArray &rd)
{
dft_log("RealArray::operator+=\n");
if (ndata != rd.ndata)
die("Size mismatch in RealArray::operator+=()\n");
for (int i=0; i < ndata; i++)
d[i] += rd.d[i];
}
/* Subtract two RealArray */
RealArray
RealArray::operator-(const RealArray &rd)
{
dft_log("RealArray::operator-\n");
if (ndata != rd.ndata)
die("Size mismatch in RealArray::operator-()\n");
RealArray rd2(*this);
for(int i = 0; i < ndata; i++)
rd2.d[i] -= rd.d[i];
return rd2;
}
/* Accumulate difference of column_bundles */
void
RealArray::operator-=(const RealArray &rd)
{
dft_log("RealArray::operator-=\n");
if (ndata != rd.ndata)
die("Size mismatch in RealArray::operator-=()\n");
for (int i=0; i < ndata; i++)
d[i] -= rd.d[i];
}
/* Scale a RealArray rd by real r */
RealArray
operator*(real r,const RealArray &rd)
{
dft_log("RealArray::operator*\n");
RealArray rd2(rd);
for(int i = 0; i < rd.ndata; i++)
rd2.d[i] *= r;
return rd2;
}
/* Scale a RealArray by real r */
RealArray
RealArray::operator*(real r)
{
dft_log("RealArray::operator*\n");
RealArray rd2(*this);
for(int i = 0; i < ndata; i++)
rd2.d[i] *= r;
return rd2;
}
/* Scale a RealArray rd by complex number c */
ComplexArray
operator*(complex c,const RealArray &rd)
{
dft_log("RealArray::operator*\n");
ComplexArray cd(rd.ndata);
for(int i =0; i < rd.ndata; i++){
cd.d[i].x = c.x*rd.d[i];
cd.d[i].y = c.y*rd.d[i];
}
return cd;
}
/* Scale a RealArray by a complex number c */
ComplexArray
RealArray::operator*(complex c)
{
ComplexArray cd(ndata);
for(int i = 0; i < ndata; i++){
cd.d[i].x = d[i]*c.x;
cd.d[i].y = d[i]*c.y;
}
return cd;
}
/* Scale a RealArray real r in place */
void
RealArray::operator*=(real r)
{
dft_log("RealArray::operator*=\n");
for (int i=0; i < ndata; i++)
d[i] *= r;
}
/////////////////////////////////////////////
// //
// MEMBER FUNCTIONS: //
// //
/////////////////////////////////////////////
// zeroes out the whole array
void
RealArray::zero_out(void)
{
for(int i=0;i < ndata; i++)
d[i] = (real)0.0;
}
/* Negates all the entries in the array */
void
RealArray::negate(void)
{
dft_log("RealArray::negate()\n");
for(int i=0;i < ndata; i++)
d[i] = -d[i];
}
/* Randomizes w/ uniform distribution in [-.5, .5]*/
void
RealArray::randomize(void)
{
//dft_log("ComplexArray::negate()\n");
for(int i=0;i < ndata; i++)
d[i] = rand()/(RAND_MAX+1.) -.5;
}
// Sum of absolute squares of all elements in rd
real abs2(const RealArray &rd)
{
real rdrd(0);
real r;
for(int i = 0; i < rd.ndata; i++){
r = rd.d[i];
rdrd += r*r;
}
return rdrd;
}
/* Take "dot-product" of two RealArray: sum the diagonals of rd1^rd2 */
real
dot(const RealArray &rd1,const RealArray &rd2)
{
dft_log("RealArray::dot()\n");
int i;
real dot12;
if (rd1.ndata != rd2.ndata)
die("rd1.ndata != rd2.ndata in dot_RealArray()\n");
dot12 = 0.0;
for (i=0; i < rd1.ndata; i++)
dot12 += rd1.d[i]*rd2.d[i];
return dot12;
}
/* Does rd2 += r * rd1 */
void
scale_accumulate(real r,const RealArray &rd1, RealArray &rd2)
{
dft_log("RealArray::scale_accumulate()\n");
for (int i=0; i < rd1.ndata; i++)
rd2.d[i] += r * rd1.d[i];
}
/* Does rd3 = r1*rd1 + r2*rd2 */
void
scaled_sum(real r1,RealArray &rd1,
real r2,RealArray &rd2,
RealArray &rd3)
{
dft_log("RealArray::scaled_sum()\n");
if (rd1.ndata != rd2.ndata || rd1.ndata != rd3.ndata)
die("Incomaptible sizes in scaled_sum(r1,rd1,r2,rd2,rd3)");
for (int i=0; i < rd1.ndata; i++)
rd3.d[i] = r1*rd1.d[i] + r2*rd2.d[i];
}
syntax highlighted by Code2HTML, v. 0.9.1