Annotation of /home/cljanss/src/SC/src/lib/math/symmetry/pointgrp.cc for ./mpqc.vmon.0018

    1                                                             //
    2                                                             // pointgrp.cc
    3                                                             //
    4                                                             // Modifications are
    5                                                             // Copyright (C) 1996 Limit Point Systems, Inc.
    6                                                             //
    7                                                             // Author: Edward Seidl <seidl@janed.com>
    8                                                             // Maintainer: LPS
    9                                                             //
   10                                                             // This file is part of the SC Toolkit.
   11                                                             //
   12                                                             // The SC Toolkit is free software; you can redistribute it and/or modify
   13                                                             // it under the terms of the GNU Library General Public License as published by
   14                                                             // the Free Software Foundation; either version 2, or (at your option)
   15                                                             // any later version.
   16                                                             //
   17                                                             // The SC Toolkit is distributed in the hope that it will be useful,
   18                                                             // but WITHOUT ANY WARRANTY; without even the implied warranty of
   19                                                             // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
   20                                                             // GNU Library General Public License for more details.
   21                                                             //
   22                                                             // You should have received a copy of the GNU Library General Public License
   23                                                             // along with the SC Toolkit; see the file COPYING.LIB.  If not, write to
   24                                                             // the Free Software Foundation, 675 Mass Ave, Cambridge, MA 02139, USA.
   25                                                             //
   26                                                             // The U.S. Government is granted a limited license as per AL 91-7.
   27                                                             //
   28                                                             
   29                                                             /* pointgrp.cc -- implementation of the point group classes
   30                                                              *
   31                                                              *      THIS SOFTWARE FITS THE DESCRIPTION IN THE U.S. COPYRIGHT ACT OF A
   32                                                              *      "UNITED STATES GOVERNMENT WORK".  IT WAS WRITTEN AS A PART OF THE
   33                                                              *      AUTHOR'S OFFICIAL DUTIES AS A GOVERNMENT EMPLOYEE.  THIS MEANS IT
   34                                                              *      CANNOT BE COPYRIGHTED.  THIS SOFTWARE IS FREELY AVAILABLE TO THE
   35                                                              *      PUBLIC FOR USE WITHOUT A COPYRIGHT NOTICE, AND THERE ARE NO
   36                                                              *      RESTRICTIONS ON ITS USE, NOW OR SUBSEQUENTLY.
   37                                                              *
   38                                                              *  Author:
   39                                                              *      E. T. Seidl
   40                                                              *      Bldg. 12A, Rm. 2033
   41                                                              *      Computer Systems Laboratory
   42                                                              *      Division of Computer Research and Technology
   43                                                              *      National Institutes of Health
   44                                                              *      Bethesda, Maryland 20892
   45                                                              *      Internet: seidl@alw.nih.gov
   46                                                              *      June, 1993
   47                                                              */
   48                                                             
   49                                                             #ifdef __GNUC__
   50                                                             #pragma implementation
   51                                                             #endif
   52                                                             
   53                                                             #include <stdlib.h>
   54                                                             #include <string.h>
   55                                                             #include <ctype.h>
   56                                                             #include <math.h>
   57                                                             
   58                                                             #include <util/misc/formio.h>
   59                                                             #include <util/state/stateio.h>
   60                                                             #include <math/symmetry/pointgrp.h>
   61                                                             
   62                                                             using namespace std;
   63                                                             
   64                                                             ////////////////////////////////////////////////////////////////////////
   65                                                             
   66                                                             static ClassDesc PointGroup_cd(
   67                                                               typeid(PointGroup),"PointGroup",2,"public SavableState",
   68                                                               create<PointGroup>, create<PointGroup>, create<PointGroup>);
   69                                                             
   70                                                             PointGroup::PointGroup()
   71                                                               : symb(0)
   72                                                             {
   73                                                               set_symbol("c1");
   74                                                               frame(0,0) = frame(1,1) = frame(2,2) = 1;
   75                                                               origin_[0] = origin_[1] = origin_[2] =0;
   76                                                             }
   77                                                             
   78                                                             PointGroup::PointGroup(const char *s)
   79                                                               : symb(0)
   80                                                             {
   81                                                               set_symbol(s);
   82                                                               frame(0,0) = frame(1,1) = frame(2,2) = 1;
   83                                                               origin_[0] = origin_[1] = origin_[2] =0;
   84                                                             }
   85                                                             
   86                                                             PointGroup::PointGroup(const char *s, SymmetryOperation& so)
   87                                                               : symb(0)
   88                                                             {
   89                                                               set_symbol(s);
   90                                                               frame = so;
   91                                                               origin_[0] = origin_[1] = origin_[2] =0;
   92                                                             }
   93                                                             
   94                                                             PointGroup::PointGroup(const char *s, SymmetryOperation& so,
   95                                                                                    const SCVector3& origin)
   96                                                               : symb(0)
   97                                                             {
   98                                                               set_symbol(s);
   99                                                               frame = so;
  100                                                               origin_ = origin;
  101                                                             }
  102                                                             
  103                                                             PointGroup::PointGroup(const Ref<KeyVal>& kv)
  104                                                               : symb(0)
  105                                                             {
  106                                                               if (kv->exists("symmetry")) {
  107                                                                 char *tmp = kv->pcharvalue("symmetry");
  108                                                                 set_symbol(tmp);
  109                                                                 delete[] tmp;
  110                                                               }
  111                                                               else
  112                                                                 set_symbol("c1");
  113                                                             
  114                                                               if (kv->exists("symmetry_frame")) {
  115                                                                 for (int i=0; i < 3; i++)
  116                                                                   for (int j=0; j < 3; j++) 
  117                                                                     frame(i,j) = kv->doublevalue("symmetry_frame",i,j);
  118                                                               } else {
  119                                                                 frame(0,0) = frame(1,1) = frame(2,2) = 1;
  120                                                               }
  121                                                             
  122                                                               if (kv->exists("origin")) {
  123                                                                 for (int i=0; i < 3; i++)
  124                                                                   origin_[i] = kv->doublevalue("origin",i);
  125                                                               } else {
  126                                                                 origin_[0] = origin_[1] = origin_[2] =0;
  127                                                               }
  128                                                             }
  129                                                             
  130                                                             PointGroup::PointGroup(StateIn& si) :
  131                                                               SavableState(si),
  132                                                               symb(0)
  133                                                             {
  134                                                               int i;
  135                                                               if (si.version(::class_desc<PointGroup>()) < 2) {
  136                                                                 ExEnv::err() << "PointGroup: checkpoint file is too old: cannot read"
  137                                                                              << endl;
  138                                                                 abort();
  139                                                               }
  140                                                               else {
  141                                                                 for (i=0; i<3; i++) si.get(origin_[i]);
  142                                                               }
  143                                                             
  144                                                               si.getstring(symb);
  145                                                               for (i=0; i < 3; i++)
  146                                                                 for (int j=0; j < 3; j++)
  147                                                                   si.get(frame(i,j));
  148                                                             }
  149                                                             
  150                                                             PointGroup::PointGroup(const PointGroup& pg)
  151                                                               : symb(0)
  152                                                             {
  153                                                               *this = pg;
  154                                                             }
  155                                                             
  156                                                             PointGroup::PointGroup(const Ref<PointGroup>& pg)
  157                                                               : symb(0)
  158                                                             {
  159                                                               *this = *pg.pointer();
  160                                                             }
  161                                                             
  162                                                             PointGroup::~PointGroup()
  163                                                             {
  164                                                               if (symb) { delete[] symb; symb=0; }
  165                                                             }
  166                                                             
  167                                                             PointGroup&
  168                                                             PointGroup::operator=(const PointGroup& pg)
  169                                                             {
  170                                                               set_symbol(pg.symb);
  171                                                               frame = pg.frame;
  172                                                               origin_ = pg.origin_;
  173                                                               return *this;
  174                                                             }
  175                                                             
  176                                                             void
  177                                                             PointGroup::set_symbol(const char *sym)
  178                                                             {
  179                                                               if (sym) {
  180                                                                 if (symb) delete[] symb;
  181                                                                 int len;
  182                                                                 symb = new char[(len=strlen(sym))+1];
  183                                                                 for (int i=0; i<len; i++) symb[i] = (char) tolower(sym[i]);
  184                                                                 symb[len] = '\0';
  185                                                               } else {
  186                                                                 set_symbol("c1");
  187                                                               }
  188                                                             }
  189                                                             
  190                                                             void
  191                                                             PointGroup::save_data_state(StateOut& so)
  192                                                             {
  193                                                               int i;
  194                                                               for (i=0; i<3; i++) so.put(origin_[i]);
  195                                                             
  196                                                               so.putstring(symb);
  197                                                             
  198                                                               for (i=0; i < 3; i++)
  199                                                                 for (int j=0; j < 3; j++)
  200                                                                   so.put(frame(i,j));
  201                                                             }
  202                                                             
  203                                                             CharacterTable
  204                                                             PointGroup::char_table() const
  205     1     0     0     0     0     0     0     0     0     0 {
  206     2     0     1     0     0     0     0     0     0     0   CharacterTable ret(symb,frame);
  207     0     0     0     1     0     0     0     0     0     0   return ret;
  208     1     0     0     0     0     0     0     0     0     0 }
  209                                                             
  210                                                             int
  211                                                             PointGroup::equiv(const Ref<PointGroup> &grp, double tol) const
  212                                                             {
  213                                                               if (strcmp(symb,grp->symb)) return 0;
  214                                                             
  215                                                               for (int i=0; i < 3; i++) {
  216                                                                 // origin isn't realy used, so don't check
  217                                                                 //if (fabs(origin_[i] - grp->origin_[i]) > tol) return 0;
  218                                                                 for (int j=0; j < 3; j++) {
  219                                                                   if (fabs(frame(i,j) - grp->frame(i,j)) > tol) return 0;
  220                                                                 }
  221                                                               }
  222                                                             
  223                                                               return 1;
  224                                                             }
  225                                                             
  226                                                             void
  227                                                             PointGroup::print(ostream &o) const
  228                                                             {
  229                                                               int i,j;
  230                                                             
  231                                                               o << node0 << indent << "symmetry = " << symb << endl;
  232                                                             
  233                                                               int unit_frame = 1;
  234                                                               int zero_origin = 1;
  235                                                               for (i=0; i<3; i++) {
  236                                                                 for (j=0; j<3; j++) {
  237                                                                   if (i==j && fabs(frame(i,j)-1.0) > 1.0e-10) unit_frame = 0;
  238                                                                   else if (i != j && fabs(frame(i,j)) > 1.0e-10) unit_frame = 0;
  239                                                                 }
  240                                                                 if (fabs(origin_[i]) > 1.0e-10) zero_origin = 0;
  241                                                               }
  242                                                             
  243                                                               if (!unit_frame) {
  244                                                                 o << node0 << indent << "symmetry_frame = [";
  245                                                                 o << incindent;
  246                                                                 for (i=0; i<3; i++) {
  247                                                                   o << node0 << endl << indent;
  248                                                                   o << node0 << "[";
  249                                                                   for (j=0; j<3; j++) {
  250                                                                     o << node0 << scprintf(" % 18.16f", frame(i,j));
  251                                                                   }
  252                                                                   o << node0 << "]";
  253                                                                 }
  254                                                                 o << node0 << "]" << endl;
  255                                                                 o << decindent;
  256                                                               }
  257                                                             
  258                                                               if (!zero_origin) {
  259                                                                 o << node0 << indent << "origin = [";
  260                                                                 for (i=0; i<3; i++) {
  261                                                                   o << node0 << scprintf(" % 18.16f", origin_[i]);
  262                                                                 }
  263                                                                 o << node0 << "]" << endl;
  264                                                               }
  265                                                             }
  266                                                             
  267                                                             /////////////////////////////////////////////////////////////////////////////
  268                                                             
  269                                                             // Local Variables:
  270                                                             // mode: c++
  271                                                             // c-file-style: "ETS"
  272                                                             // End:
  273