FOX/ObjCryst++  1.10.X (development)
1 /* ObjCryst++ Object-Oriented Crystallographic Library
2  (c) 2000-2002 Vincent Favre-Nicolin
3  2000-2001 University of Geneva (Switzerland)
5  This program is free software; you can redistribute it and/or modify
6  it under the terms of the GNU General Public License as published by
7  the Free Software Foundation; version 2 of the License.
9  This program is distributed in the hope that it will be useful,
10  but WITHOUT ANY WARRANTY; without even the implied warranty of
12  GNU General Public License for more details.
14  You should have received a copy of the GNU General Public License
15  along with this program; if not, write to the Free Software
16  Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
17 */
18 /* test.cpp
19 * source file for test functions (speed, etc...)
20 *
21 */
22 #include <stdlib.h>
23 #include <list>
24 #include "ObjCryst/ObjCryst/test.h"
25 #include "ObjCryst/ObjCryst/Crystal.h"
26 #include "ObjCryst/ObjCryst/Atom.h"
27 #include "ObjCryst/ObjCryst/DiffractionDataSingleCrystal.h"
28 #include "ObjCryst/ObjCryst/PowderPattern.h"
29 #include "ObjCryst/RefinableObj/GlobalOptimObj.h"
30 #include "ObjCryst/Quirks/VFNStreamFormat.h"
32 namespace ObjCryst
33 {
35 SpeedTestReport SpeedTest(const unsigned int nbAtom, const int nbAtomType,const string spacegroup,
36  const RadiationType radiation, const unsigned long nbReflections,
37  const unsigned int dataType,const REAL time)
38 {
39  Crystal cryst(9,11,15,1.2,1.3,1.7,spacegroup);
40  for(int i=0;i<nbAtomType;++i)
41  {
42  cryst.AddScatteringPower(new ScatteringPowerAtom("O","O",1.5));
43  }
44  for(unsigned int i = 0; i < nbAtom; ++i)
45  {
46  cryst.AddScatterer(new Atom(.0,.0,.0,"O",
47  &(cryst.GetScatteringPowerRegistry().GetObj(i%nbAtomType)),
48  1.));
49  }
50  cryst.SetUseDynPopCorr(false);
52  RefinableObj *pData = NULL;
53  switch(dataType)
54  {
55  case 0:
56  {
58  pDataTmp->SetWavelength(0.25);
59  pDataTmp->SetRadiationType(radiation);
60  pDataTmp->SetMaxSinThetaOvLambda(100.);
61  pDataTmp->SetCrystal(cryst);
62  float maxtheta=0.1;
63  for(;;)
64  {
65  pDataTmp->GenHKLFullSpace(maxtheta, true);
66  if(pDataTmp->GetNbRefl()>(long)nbReflections) break;
67  maxtheta*=1.5;
68  if(maxtheta>=M_PI/2.) break;
69  }
70  CrystVector_REAL hh; hh=pDataTmp->GetH();hh.resizeAndPreserve(nbReflections);hh+=0.0001;
71  CrystVector_REAL kk; kk=pDataTmp->GetK();kk.resizeAndPreserve(nbReflections);kk+=0.0001;
72  CrystVector_REAL ll; ll=pDataTmp->GetL();ll.resizeAndPreserve(nbReflections);ll+=0.0001;
74  CrystVector_long h(nbReflections); h=hh;
75  CrystVector_long k(nbReflections); k=kk;
76  CrystVector_long l(nbReflections); l=ll;
78  CrystVector_REAL iobs(nbReflections);
79  for(unsigned int i=0;i<nbReflections;++i) iobs(i)=(REAL)rand();
80  CrystVector_REAL sigma(nbReflections);sigma=1;
82  pDataTmp->SetHklIobs (h, k, l, iobs, sigma);
83  pDataTmp->SetWeightToInvSigma2();
85  pData=pDataTmp;
86  break;
87  }
88  case 1:
89  {
90  PowderPattern *pDataTmp=new PowderPattern;
91  pDataTmp->SetWavelength(0.25);
92  pDataTmp->SetRadiationType(radiation);
93  pDataTmp->SetPowderPatternPar(0.001,.001,3140);
94  CrystVector_REAL iobs(3140);
95  for(unsigned int i=0;i<3140;++i) iobs(i)=(REAL)rand()+1.;
96  pDataTmp->SetPowderPatternObs(iobs);
97  pDataTmp->SetMaxSinThetaOvLambda(100.);
100  backgdData->SetName("PbSo4-background");
101  {
102  CrystVector_REAL tth(2),backgd(2);
103  tth(0)=0.;tth(1)=3.14;
104  backgd(0)=1.;backgd(1)=9.;
105  backgdData->SetInterpPoints(tth,backgd);
106  }
107  pDataTmp->AddPowderPatternComponent(*backgdData);
110  diffData->SetCrystal(cryst);
111  pDataTmp->AddPowderPatternComponent(*diffData);
112  diffData->SetName("Crystal phase");
113  diffData->SetReflectionProfilePar(PROFILE_PSEUDO_VOIGT,
114  .03*DEG2RAD*DEG2RAD,0.,0.,0.3,0);
115  {
116  float maxtheta=0.1;
117  for(;;)
118  {
119  diffData->ScatteringData::GenHKLFullSpace(maxtheta, true);
120  if(diffData->GetNbRefl()>(long)nbReflections) break;
121  maxtheta*=1.5;
122  if(maxtheta>=M_PI/2.) break;
123  }
124  CrystVector_REAL hh; hh=diffData->GetH();hh.resizeAndPreserve(nbReflections);
125  CrystVector_REAL kk; kk=diffData->GetK();kk.resizeAndPreserve(nbReflections);
126  CrystVector_REAL ll; ll=diffData->GetL();ll.resizeAndPreserve(nbReflections);
128  diffData->SetHKL (hh, kk, ll);
129  }
130  pData=pDataTmp;
131  break;
132  }
133  }
135  //Create the global optimization object
136  MonteCarloObj *pGlobalOptObj=new MonteCarloObj;
137  pGlobalOptObj->AddRefinableObj(*pData);
138  pGlobalOptObj->AddRefinableObj(cryst);
140  //Refine only positionnal parameters
141  pGlobalOptObj->FixAllPar();
142  pGlobalOptObj->SetParIsFixed(gpRefParTypeScattTransl,false);
143  pGlobalOptObj->SetParIsFixed(gpRefParTypeScattOrient,false);
145  //Don't cheat ;-)
146  pGlobalOptObj->RandomizeStartingConfig();
148  //Annealing parameters (schedule, Tmax, Tmin, displacement schedule,
149  pGlobalOptObj->SetAlgorithmParallTempering(ANNEALING_SMART,1e8,1e-8,
152  //Global Optimization
153  //The real job-first test
154  long nbTrial=50000000;
155  pGlobalOptObj->Optimize(nbTrial,true,0,time);
158  SpeedTestReport report;
159  report.mNbAtom=nbAtom;
160  report.mNbAtomType=nbAtomType;
161  report.mSpacegroup=spacegroup;
162  report.mRadiation=radiation;
163  report.mNbReflections=nbReflections;
164  report.mDataType=dataType;
165  report.mBogoMRAPS=(REAL)nbAtom*cryst.GetSpaceGroup().GetNbSymmetrics()*(REAL)nbReflections
166  *(50000000-nbTrial)/pGlobalOptObj->GetLastOptimElapsedTime()/1e6;
167  report.mBogoMRAPS_reduced=(REAL)nbAtom*cryst.GetSpaceGroup().GetNbSymmetrics(true,true)
168  *(REAL)nbReflections
169  *(50000000-nbTrial)/pGlobalOptObj->GetLastOptimElapsedTime()/1e6;
170  report.mBogoSPS=(50000000-nbTrial)/pGlobalOptObj->GetLastOptimElapsedTime();
171  delete pGlobalOptObj;
172  delete pData;
173  return report;
174 }
175 }
const CrystVector_REAL & GetH() const
Return the 1D array of H coordinates for all reflections.
int GetNbSymmetrics(const bool noCenter=false, const bool noTransl=false) const
Return the number of equivalent positions in the spacegroup, ie the multilicity of the general positi...
Definition: SpaceGroup.cpp:419
virtual void RandomizeStartingConfig()
Randomize starting configuration.
virtual void Optimize(long &nbSteps, const bool silent=false, const REAL finalcost=0, const REAL maxTime=-1)
Launch optimization (a single run) for N steps.
unsigned int mDataType
dataType: 0= single crystal, 1= powder pattern (1 background + 1 crystal phase)
Definition: test.h:46
void SetReflectionProfilePar(const ReflectionProfileType prof, const REAL fwhmCagliotiW, const REAL fwhmCagliotiU=0, const REAL fwhmCagliotiV=0, const REAL eta0=0.5, const REAL eta1=0.)
Set reflection profile parameters.
Million of Reflections-Atoms computed Per Second (considering all atoms in the unit cell) ...
Definition: test.h:48
void SetUseDynPopCorr(const int use)
Set the use of dynamical population correction (Crystal::mUseDynPopCorr).
Definition: Crystal.cpp:791
void AddRefinableObj(RefinableObj &)
Add a refined object. All sub-objects are also added.
RadiationType mRadiation
The type of radiation used.
Definition: test.h:42
Phase to compute a background contribution to a powder pattern using an interpolation.
void SetParIsFixed(const string &parName, const bool fix)
Fix one parameter.
virtual void SetCrystal(Crystal &crystal)
Set the crystal for this experiment.
Class to compute the contribution to a powder pattern from a crystalline phase.
REAL GetLastOptimElapsedTime() const
Get the elapsed time (in seconds) during the last optimization.
void SetPowderPatternObs(const CrystVector_REAL &obs)
Set observed powder pattern from vector array.
void FixAllPar()
Fix all parameters.
virtual void SetMaxSinThetaOvLambda(const REAL max)
Set the maximum value for sin(theta)/lambda.
unsigned int mNbAtom
Total number of unique atoms in the test structure.
Definition: test.h:36
const CrystVector_REAL & GetL() const
Return the 1D array of L coordinates for all reflections.
Generic Refinable Object.
Definition: RefinableObj.h:752
const SpaceGroup & GetSpaceGroup() const
Access to the SpaceGroup object.
Definition: UnitCell.cpp:320
unsigned long mNbReflections
The total number of reflections used for the tests.
Definition: test.h:44
void AddScatterer(Scatterer *scatt)
Add a scatterer to the crystal.
Definition: Crystal.cpp:174
REAL mBogoMRAPS_reduced
Million of Reflections-Atoms computed Per Second (considering all atoms in the unit cell...
Definition: test.h:51
virtual void GenHKLFullSpace(const REAL maxTheta, const bool unique=false)
Generate a list of h,k,l to describe a full reciprocal space, up to a given maximum theta value...
Base object for Monte-Carlo Global Optimization methods.
The basic atom scatterer, in a crystal.
Definition: Atom.h:57
virtual void SetRadiationType(const RadiationType radiation)
Set : neutron or x-ray experiment ? Wavelength ?
void SetAlgorithmParallTempering(const AnnealingSchedule scheduleTemp, const REAL tMax, const REAL tMin, const AnnealingSchedule scheduleMutation=ANNEALING_CONSTANT, const REAL mutMax=16., const REAL mutMin=.125)
Set the refinement method to Parallel Tempering.
DiffractionData object for Single Crystal analysis.
void AddPowderPatternComponent(PowderPatternComponent &)
Add a component (phase, backround) to this pattern.
string mSpacegroup
The symbol for the spacegroup.
Definition: test.h:40
virtual void SetMaxSinThetaOvLambda(const REAL max)
Set the maximum value for sin(theta)/lambda.
void SetHklIobs(const CrystVector_long &h, const CrystVector_long &k, const CrystVector_long &l, const CrystVector_REAL &iObs, const CrystVector_REAL &sigma)
input H,K,L, Iobs and Sigma
Powder pattern class, with an observed pattern and several calculated components to modelize the patt...
void SetPowderPatternPar(const REAL min, const REAL step, unsigned long nbPoint)
the powder pattern angular range & resolution parameter.
void SetRadiationType(const RadiationType radiation)
Set the radiation type.
virtual void SetCrystal(Crystal &crystal)
Set the crystal for this experiment.
Number of Structures evaluated Per Second.
Definition: test.h:53
ObjRegistry< ScatteringPower > & GetScatteringPowerRegistry()
Get the registry of ScatteringPower included in this Crystal.
Definition: Crystal.cpp:222
Structure to hold the results of a speedtest (see ObjCryst::SpeedTest())
Definition: test.h:33
int mNbAtomType
Total number of atom types in the test structure.
Definition: test.h:38
SpeedTestReport SpeedTest(const unsigned int nbAtom, const int nbAtomType, const string spacegroup, const RadiationType radiation, const unsigned long nbReflections, const unsigned int dataType, const REAL time)
Definition: test.cpp:35
virtual void SetWeightToInvSigma2(const REAL minRelatSigma=1e-4, const REAL minIobsSigmaRatio=0)
Set the weight for all observed intensities to 1/sigma^2.
long GetNbRefl() const
Return the number of reflections in this experiment.
The namespace which includes all objects (crystallographic and algorithmic) in ObjCryst++.
Definition: Atom.cpp:47
The Scattering Power for an Atom.
void SetWavelength(const REAL lambda)
Set the wavelength of the experiment (in Angstroems).
virtual void SetHKL(const CrystVector_REAL &h, const CrystVector_REAL &k, const CrystVector_REAL &l)
input H,K,L
Type of radiation used.
Definition: General.h:94
Crystal class: Unit cell, spacegroup, scatterers.
Definition: Crystal.h:97
const CrystVector_REAL & GetK() const
Return the 1D array of K coordinates for all reflections.
void AddScatteringPower(ScatteringPower *scattPow)
Add a ScatteringPower for this Crystal.
Definition: Crystal.cpp:227
virtual void SetName(const string &name)
Name of the object.
void SetWavelength(const REAL)
Set the (monochromatic) wavelength of the beam.