testBasis.C 2.42 KB
Newer Older
Francois Gygi committed
1 2
////////////////////////////////////////////////////////////////////////////////
//
Francois Gygi committed
3 4 5 6
// Copyright (c) 2008 The Regents of the University of California
//
// This file is part of Qbox
//
Francois Gygi committed
7 8
// Qbox is distributed under the terms of the GNU General Public License
// as published by the Free Software Foundation, either version 2 of
Francois Gygi committed
9 10 11 12 13 14
// the License, or (at your option) any later version.
// See the file COPYING in the root directory of this distribution
// or <http://www.gnu.org/licenses/>.
//
////////////////////////////////////////////////////////////////////////////////
//
Francois Gygi committed
15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36
// testBasis.C
//
////////////////////////////////////////////////////////////////////////////////

#include "Basis.h"
#include <iostream>
#include <new>
#include <cstdlib>
#include <cassert>
using namespace std;

#ifdef USE_MPI
#include <mpi.h>
#endif

int main(int argc, char **argv)
{
  // use: testBasis a0x a0y a0z a1x a1y a1z a2x a2y a2z ecut kx ky kz npr npc
#if USE_MPI
  MPI_Init(&argc,&argv);
#endif
  {
Francois Gygi committed
37 38 39 40 41 42 43
    if ( argc !=16 )
    {
      cout <<
      " use: testBasis a0x a0y a0z a1x a1y a1z a2x a2y a2z ecut kx ky kz npr npc"
      << endl;
      return 1;
    }
Francois Gygi committed
44 45 46
    const D3vector a0(atof(argv[1]),atof(argv[2]),atof(argv[3]));
    const D3vector a1(atof(argv[4]),atof(argv[5]),atof(argv[6]));
    const D3vector a2(atof(argv[7]),atof(argv[8]),atof(argv[9]));
47

Francois Gygi committed
48 49 50 51
    double ecut = atof(argv[10]);
    D3vector kpoint(atof(argv[11]),atof(argv[12]),atof(argv[13]));
    int npr = atoi(argv[14]);
    int npc = atoi(argv[15]);
52

Francois Gygi committed
53
    UnitCell cell(a0,a1,a2);
54

55 56 57 58 59 60 61 62 63 64 65 66 67 68
    // create cartesian communicator
    int ndims=2;
    int dims[2] = {npr, npc};
    int periods[2] = { 0, 0};
    int reorder = 0;
    MPI_Comm cartcomm;
    MPI_Cart_create(MPI_COMM_WORLD,ndims,dims,periods,reorder,&cartcomm);

    // partition the cartesian communicator in columns
    MPI_Comm colcomm;
    int remain_dims[2] = { 1, 0 };
    MPI_Cart_sub(cartcomm,remain_dims,&colcomm);

    Basis basis(colcomm,kpoint);
Francois Gygi committed
69 70 71 72 73 74 75 76 77
    try
    {
      basis.resize(cell,cell,ecut);
    }
    catch ( bad_alloc )
    {
      cout << " bad_alloc caught in Basis::resize" << endl;
      throw;
    }
78

79 80 81 82
    int npes, mype;
    MPI_Comm_size(colcomm,&npes);
    MPI_Comm_rank(colcomm,&mype);
    for (int ipe = 0; ipe < npes; ipe++ )
83
    {
84 85
      MPI_Barrier(colcomm);
      if ( ipe == mype )
86 87 88 89 90
      {
        cout << basis;
        cout << endl;
      }
    }
91

Francois Gygi committed
92 93 94 95 96 97 98
    //Basis b2(basis);
    //cout << b2;
  }
#if USE_MPI
  MPI_Finalize();
#endif
}