source: src/LinearAlgebra/unittests/MatrixUnitTest.cpp@ 311da7b

Action_Thermostats Add_AtomRandomPerturbation Add_FitFragmentPartialChargesAction Add_RotateAroundBondAction Add_SelectAtomByNameAction Added_ParseSaveFragmentResults AddingActions_SaveParseParticleParameters Adding_Graph_to_ChangeBondActions Adding_MD_integration_tests Adding_ParticleName_to_Atom Adding_StructOpt_integration_tests AtomFragments Automaking_mpqc_open AutomationFragmentation_failures Candidate_v1.5.4 Candidate_v1.6.0 Candidate_v1.6.1 ChangeBugEmailaddress ChangingTestPorts ChemicalSpaceEvaluator CombiningParticlePotentialParsing Combining_Subpackages Debian_Package_split Debian_package_split_molecuildergui_only Disabling_MemDebug Docu_Python_wait EmpiricalPotential_contain_HomologyGraph EmpiricalPotential_contain_HomologyGraph_documentation Enable_parallel_make_install Enhance_userguide Enhanced_StructuralOptimization Enhanced_StructuralOptimization_continued Example_ManyWaysToTranslateAtom Exclude_Hydrogens_annealWithBondGraph FitPartialCharges_GlobalError Fix_BoundInBox_CenterInBox_MoleculeActions Fix_ChargeSampling_PBC Fix_ChronosMutex Fix_FitPartialCharges Fix_FitPotential_needs_atomicnumbers Fix_ForceAnnealing Fix_IndependentFragmentGrids Fix_ParseParticles Fix_ParseParticles_split_forward_backward_Actions Fix_PopActions Fix_QtFragmentList_sorted_selection Fix_Restrictedkeyset_FragmentMolecule Fix_StatusMsg Fix_StepWorldTime_single_argument Fix_Verbose_Codepatterns Fix_fitting_potentials Fixes ForceAnnealing_goodresults ForceAnnealing_oldresults ForceAnnealing_tocheck ForceAnnealing_with_BondGraph ForceAnnealing_with_BondGraph_continued ForceAnnealing_with_BondGraph_continued_betteresults ForceAnnealing_with_BondGraph_contraction-expansion FragmentAction_writes_AtomFragments FragmentMolecule_checks_bonddegrees GeometryObjects Gui_Fixes Gui_displays_atomic_force_velocity ImplicitCharges IndependentFragmentGrids IndependentFragmentGrids_IndividualZeroInstances IndependentFragmentGrids_IntegrationTest IndependentFragmentGrids_Sole_NN_Calculation JobMarket_RobustOnKillsSegFaults JobMarket_StableWorkerPool JobMarket_unresolvable_hostname_fix MoreRobust_FragmentAutomation ODR_violation_mpqc_open PartialCharges_OrthogonalSummation PdbParser_setsAtomName PythonUI_with_named_parameters QtGui_reactivate_TimeChanged_changes Recreated_GuiChecks Rewrite_FitPartialCharges RotateToPrincipalAxisSystem_UndoRedo SaturateAtoms_findBestMatching SaturateAtoms_singleDegree StoppableMakroAction Subpackage_CodePatterns Subpackage_JobMarket Subpackage_LinearAlgebra Subpackage_levmar Subpackage_mpqc_open Subpackage_vmg Switchable_LogView ThirdParty_MPQC_rebuilt_buildsystem TrajectoryDependenant_MaxOrder TremoloParser_IncreasedPrecision TremoloParser_MultipleTimesteps TremoloParser_setsAtomName Ubuntu_1604_changes stable
Last change on this file since 311da7b was fff54f, checked in by Frederik Heber <heber@…>, 14 years ago

Moved all Matrix...Unittest to LinearAlgebra/unittests/Matrix...UnitTest.

  • Property mode set to 100644
File size: 10.4 KB
RevLine 
[bcf653]1/*
2 * Project: MoleCuilder
3 * Description: creates and alters molecular systems
4 * Copyright (C) 2010 University of Bonn. All rights reserved.
5 * Please see the LICENSE file or "Copyright notice" in builder.cpp for details.
6 */
7
[5630bd]8/*
[fff54f]9 * MatrixUnitTest.cpp
[5630bd]10 *
11 * Created on: Jul 7, 2010
12 * Author: crueger
13 */
14
[bf3817]15// include config.h
16#ifdef HAVE_CONFIG_H
17#include <config.h>
18#endif
19
[5630bd]20#include <cppunit/CompilerOutputter.h>
21#include <cppunit/extensions/TestFactoryRegistry.h>
22#include <cppunit/ui/text/TestRunner.h>
23
[31fb1d]24#include <cmath>
25
[fff54f]26#include "MatrixUnitTest.hpp"
[57f243]27#include "LinearAlgebra/Vector.hpp"
[cca9ef]28#include "LinearAlgebra/RealSpaceMatrix.hpp"
[5630bd]29#include "Exceptions/NotInvertibleException.hpp"
30
31#ifdef HAVE_TESTRUNNER
32#include "UnitTestMain.hpp"
33#endif /*HAVE_TESTRUNNER*/
34
35// Registers the fixture into the 'registry'
36CPPUNIT_TEST_SUITE_REGISTRATION( MatrixUnittest );
37
38void MatrixUnittest::setUp(){
[cca9ef]39 zero = new RealSpaceMatrix();
[3bc926]40 for(int i =NDIM;i--;) {
41 for(int j =NDIM;j--;) {
42 zero->at(i,j)=0.;
43 }
44 }
[cca9ef]45 one = new RealSpaceMatrix();
[5630bd]46 for(int i =NDIM;i--;){
47 one->at(i,i)=1.;
48 }
[cca9ef]49 full=new RealSpaceMatrix();
[5630bd]50 for(int i=NDIM;i--;){
51 for(int j=NDIM;j--;){
52 full->at(i,j)=1.;
53 }
54 }
[cca9ef]55 diagonal = new RealSpaceMatrix();
[5630bd]56 for(int i=NDIM;i--;){
57 diagonal->at(i,i)=i+1.;
58 }
[cca9ef]59 perm1 = new RealSpaceMatrix();
[407782]60 perm1->column(0) = unitVec[0];
61 perm1->column(1) = unitVec[2];
62 perm1->column(2) = unitVec[1];
[5630bd]63
64
[cca9ef]65 perm2 = new RealSpaceMatrix();
[407782]66 perm2->column(0) = unitVec[1];
67 perm2->column(1) = unitVec[0];
68 perm2->column(2) = unitVec[2];
[5630bd]69
[cca9ef]70 perm3 = new RealSpaceMatrix();
[407782]71 perm3->column(0) = unitVec[1];
72 perm3->column(1) = unitVec[2];
73 perm3->column(2) = unitVec[0];
[5630bd]74
[cca9ef]75 perm4 = new RealSpaceMatrix();
[407782]76 perm4->column(0) = unitVec[2];
77 perm4->column(1) = unitVec[1];
78 perm4->column(2) = unitVec[0];
[5630bd]79
[cca9ef]80 perm5 = new RealSpaceMatrix();
[407782]81 perm5->column(0) = unitVec[2];
82 perm5->column(1) = unitVec[0];
83 perm5->column(2) = unitVec[1];
[5630bd]84
85}
86void MatrixUnittest::tearDown(){
87 delete zero;
88 delete one;
89 delete full;
90 delete diagonal;
91 delete perm1;
92 delete perm2;
93 delete perm3;
94 delete perm4;
95 delete perm5;
96}
97
98void MatrixUnittest::AccessTest(){
[cca9ef]99 RealSpaceMatrix mat;
[5630bd]100 for(int i=NDIM;i--;){
101 for(int j=NDIM;j--;){
102 CPPUNIT_ASSERT_EQUAL(mat.at(i,j),0.);
103 }
104 }
105 int k=1;
106 for(int i=NDIM;i--;){
107 for(int j=NDIM;j--;){
108 mat.at(i,j)=k++;
109 }
110 }
111 k=1;
112 for(int i=NDIM;i--;){
113 for(int j=NDIM;j--;){
114 CPPUNIT_ASSERT_EQUAL(mat.at(i,j),(double)k);
115 ++k;
116 }
117 }
118}
119
120void MatrixUnittest::VectorTest(){
[cca9ef]121 RealSpaceMatrix mat;
[5630bd]122 for(int i=NDIM;i--;){
123 CPPUNIT_ASSERT_EQUAL(mat.row(i),zeroVec);
124 CPPUNIT_ASSERT_EQUAL(mat.column(i),zeroVec);
125 }
126 CPPUNIT_ASSERT_EQUAL(mat.diagonal(),zeroVec);
127
[9eb7580]128 mat.setIdentity();
[407782]129 CPPUNIT_ASSERT_EQUAL(mat.row(0),unitVec[0]);
130 CPPUNIT_ASSERT_EQUAL(mat.row(1),unitVec[1]);
131 CPPUNIT_ASSERT_EQUAL(mat.row(2),unitVec[2]);
132 CPPUNIT_ASSERT_EQUAL(mat.column(0),unitVec[0]);
133 CPPUNIT_ASSERT_EQUAL(mat.column(1),unitVec[1]);
134 CPPUNIT_ASSERT_EQUAL(mat.column(2),unitVec[2]);
[5630bd]135
136 Vector t1=Vector(1.,1.,1.);
137 Vector t2=Vector(2.,2.,2.);
138 Vector t3=Vector(3.,3.,3.);
139 Vector t4=Vector(1.,2.,3.);
140
141 mat.row(0)=t1;
142 mat.row(1)=t2;
143 mat.row(2)=t3;
144 CPPUNIT_ASSERT_EQUAL(mat.row(0),t1);
145 CPPUNIT_ASSERT_EQUAL(mat.row(1),t2);
146 CPPUNIT_ASSERT_EQUAL(mat.row(2),t3);
147 CPPUNIT_ASSERT_EQUAL(mat.column(0),t4);
148 CPPUNIT_ASSERT_EQUAL(mat.column(1),t4);
149 CPPUNIT_ASSERT_EQUAL(mat.column(2),t4);
150 CPPUNIT_ASSERT_EQUAL(mat.diagonal(),t4);
151 for(int i=NDIM;i--;){
152 for(int j=NDIM;j--;){
153 CPPUNIT_ASSERT_EQUAL(mat.at(i,j),i+1.);
154 }
155 }
156
157 mat.column(0)=t1;
158 mat.column(1)=t2;
159 mat.column(2)=t3;
160 CPPUNIT_ASSERT_EQUAL(mat.column(0),t1);
161 CPPUNIT_ASSERT_EQUAL(mat.column(1),t2);
162 CPPUNIT_ASSERT_EQUAL(mat.column(2),t3);
163 CPPUNIT_ASSERT_EQUAL(mat.row(0),t4);
164 CPPUNIT_ASSERT_EQUAL(mat.row(1),t4);
165 CPPUNIT_ASSERT_EQUAL(mat.row(2),t4);
166 CPPUNIT_ASSERT_EQUAL(mat.diagonal(),t4);
167 for(int i=NDIM;i--;){
168 for(int j=NDIM;j--;){
169 CPPUNIT_ASSERT_EQUAL(mat.at(i,j),j+1.);
170 }
171 }
172}
173
[31fb1d]174void MatrixUnittest::TransposeTest(){
[cca9ef]175 RealSpaceMatrix res;
[31fb1d]176
177 // transpose of unit is unit
[9eb7580]178 res.setIdentity();
[3bc926]179 res.transpose();
[31fb1d]180 CPPUNIT_ASSERT_EQUAL(res,*one);
181
182 // transpose of transpose is same matrix
[9eb7580]183 res.setZero();
[31fb1d]184 res.set(2,2, 1.);
185 CPPUNIT_ASSERT_EQUAL(res.transpose().transpose(),res);
186}
187
[5630bd]188void MatrixUnittest::OperationTest(){
[cca9ef]189 RealSpaceMatrix res;
[5630bd]190
191 res =(*zero) *(*zero);
[3bc926]192 std::cout << *zero << " times " << *zero << " is " << res << std::endl;
[5630bd]193 CPPUNIT_ASSERT_EQUAL(res,*zero);
194 res =(*zero) *(*one);
195 CPPUNIT_ASSERT_EQUAL(res,*zero);
196 res =(*zero) *(*full);
197 CPPUNIT_ASSERT_EQUAL(res,*zero);
198 res =(*zero) *(*diagonal);
199 CPPUNIT_ASSERT_EQUAL(res,*zero);
200 res =(*zero) *(*perm1);
201 CPPUNIT_ASSERT_EQUAL(res,*zero);
202 res =(*zero) *(*perm2);
203 CPPUNIT_ASSERT_EQUAL(res,*zero);
204 res =(*zero) *(*perm3);
205 CPPUNIT_ASSERT_EQUAL(res,*zero);
206 res =(*zero) *(*perm4);
207 CPPUNIT_ASSERT_EQUAL(res,*zero);
208 res =(*zero) *(*perm5);
209 CPPUNIT_ASSERT_EQUAL(res,*zero);
210
211 res =(*one)*(*one);
212 CPPUNIT_ASSERT_EQUAL(res,*one);
213 res =(*one)*(*full);
214 CPPUNIT_ASSERT_EQUAL(res,*full);
215 res =(*one)*(*diagonal);
216 CPPUNIT_ASSERT_EQUAL(res,*diagonal);
217 res =(*one)*(*perm1);
218 CPPUNIT_ASSERT_EQUAL(res,*perm1);
219 res =(*one)*(*perm2);
220 CPPUNIT_ASSERT_EQUAL(res,*perm2);
221 res =(*one)*(*perm3);
222 CPPUNIT_ASSERT_EQUAL(res,*perm3);
223 res =(*one)*(*perm4);
224 CPPUNIT_ASSERT_EQUAL(res,*perm4);
225 res =(*one)*(*perm5);
226 CPPUNIT_ASSERT_EQUAL(res,*perm5);
227
228 res = (*full)*(*perm1);
229 CPPUNIT_ASSERT_EQUAL(res,*full);
230 res = (*full)*(*perm2);
231 CPPUNIT_ASSERT_EQUAL(res,*full);
232 res = (*full)*(*perm3);
233 CPPUNIT_ASSERT_EQUAL(res,*full);
234 res = (*full)*(*perm4);
235 CPPUNIT_ASSERT_EQUAL(res,*full);
236 res = (*full)*(*perm5);
237 CPPUNIT_ASSERT_EQUAL(res,*full);
238
239 res = (*diagonal)*(*perm1);
[407782]240 CPPUNIT_ASSERT_EQUAL(res.column(0),unitVec[0]);
241 CPPUNIT_ASSERT_EQUAL(res.column(1),3*unitVec[2]);
242 CPPUNIT_ASSERT_EQUAL(res.column(2),2*unitVec[1]);
[5630bd]243 res = (*diagonal)*(*perm2);
[407782]244 CPPUNIT_ASSERT_EQUAL(res.column(0),2*unitVec[1]);
245 CPPUNIT_ASSERT_EQUAL(res.column(1),unitVec[0]);
246 CPPUNIT_ASSERT_EQUAL(res.column(2),3*unitVec[2]);
[5630bd]247 res = (*diagonal)*(*perm3);
[407782]248 CPPUNIT_ASSERT_EQUAL(res.column(0),2*unitVec[1]);
249 CPPUNIT_ASSERT_EQUAL(res.column(1),3*unitVec[2]);
250 CPPUNIT_ASSERT_EQUAL(res.column(2),unitVec[0]);
[5630bd]251 res = (*diagonal)*(*perm4);
[407782]252 CPPUNIT_ASSERT_EQUAL(res.column(0),3*unitVec[2]);
253 CPPUNIT_ASSERT_EQUAL(res.column(1),2*unitVec[1]);
254 CPPUNIT_ASSERT_EQUAL(res.column(2),unitVec[0]);
[5630bd]255 res = (*diagonal)*(*perm5);
[407782]256 CPPUNIT_ASSERT_EQUAL(res.column(0),3*unitVec[2]);
257 CPPUNIT_ASSERT_EQUAL(res.column(1),unitVec[0]);
258 CPPUNIT_ASSERT_EQUAL(res.column(2),2*unitVec[1]);
[5630bd]259}
260
[31fb1d]261void MatrixUnittest::RotationTest(){
[cca9ef]262 RealSpaceMatrix res;
263 RealSpaceMatrix inverse;
[31fb1d]264
265 // zero rotation angles yields unity matrix
[9eb7580]266 res.setRotation(0,0,0);
[31fb1d]267 CPPUNIT_ASSERT_EQUAL(*one, res);
268
269 // arbitrary rotation matrix has det = 1
[9eb7580]270 res.setRotation(M_PI/3.,1.,M_PI/7.);
[31fb1d]271 CPPUNIT_ASSERT(fabs(fabs(res.determinant()) -1.) < MYEPSILON);
272
273 // inverse is rotation matrix with negative angles
[9eb7580]274 res.setRotation(M_PI/3.,0.,0.);
275 inverse.setRotation(-M_PI/3.,0.,0.);
[31fb1d]276 CPPUNIT_ASSERT_EQUAL(*one, res * inverse);
277
278 // ... or transposed
[9eb7580]279 res.setRotation(M_PI/3.,0.,0.);
[cca9ef]280 CPPUNIT_ASSERT_EQUAL(inverse, ((const RealSpaceMatrix) res).transpose());
[31fb1d]281}
[5630bd]282
283void MatrixUnittest::InvertTest(){
284 CPPUNIT_ASSERT_THROW(zero->invert(),NotInvertibleException);
285 CPPUNIT_ASSERT_THROW(full->invert(),NotInvertibleException);
286
[cca9ef]287 RealSpaceMatrix res;
[5630bd]288 res = (*one)*one->invert();
289 CPPUNIT_ASSERT_EQUAL(res,*one);
290 res = (*diagonal)*diagonal->invert();
291 CPPUNIT_ASSERT_EQUAL(res,*one);
292 res = (*perm1)*perm1->invert();
293 CPPUNIT_ASSERT_EQUAL(res,*one);
294 res = (*perm2)*perm2->invert();
295 CPPUNIT_ASSERT_EQUAL(res,*one);
296 res = (*perm3)*perm3->invert();
297 CPPUNIT_ASSERT_EQUAL(res,*one);
298 res = (*perm4)*perm4->invert();
299 CPPUNIT_ASSERT_EQUAL(res,*one);
300 res = (*perm5)*perm5->invert();
301 CPPUNIT_ASSERT_EQUAL(res,*one);
302}
303
304
305void MatrixUnittest::DeterminantTest(){
306 CPPUNIT_ASSERT_EQUAL(zero->determinant(),0.);
307 CPPUNIT_ASSERT_EQUAL(one->determinant(),1.);
308 CPPUNIT_ASSERT_EQUAL(diagonal->determinant(),6.);
309 CPPUNIT_ASSERT_EQUAL(full->determinant(),0.);
310 CPPUNIT_ASSERT_EQUAL(perm1->determinant(),-1.);
311 CPPUNIT_ASSERT_EQUAL(perm2->determinant(),-1.);
312 CPPUNIT_ASSERT_EQUAL(perm3->determinant(),1.);
313 CPPUNIT_ASSERT_EQUAL(perm4->determinant(),-1.);
314 CPPUNIT_ASSERT_EQUAL(perm5->determinant(),1.);
315}
316
317void MatrixUnittest::VecMultTest(){
[407782]318 CPPUNIT_ASSERT_EQUAL((*zero)*unitVec[0],zeroVec);
319 CPPUNIT_ASSERT_EQUAL((*zero)*unitVec[1],zeroVec);
320 CPPUNIT_ASSERT_EQUAL((*zero)*unitVec[2],zeroVec);
[5630bd]321 CPPUNIT_ASSERT_EQUAL((*zero)*zeroVec,zeroVec);
322
[407782]323 CPPUNIT_ASSERT_EQUAL((*one)*unitVec[0],unitVec[0]);
324 CPPUNIT_ASSERT_EQUAL((*one)*unitVec[1],unitVec[1]);
325 CPPUNIT_ASSERT_EQUAL((*one)*unitVec[2],unitVec[2]);
[5630bd]326 CPPUNIT_ASSERT_EQUAL((*one)*zeroVec,zeroVec);
327
[407782]328 CPPUNIT_ASSERT_EQUAL((*diagonal)*unitVec[0],unitVec[0]);
329 CPPUNIT_ASSERT_EQUAL((*diagonal)*unitVec[1],2*unitVec[1]);
330 CPPUNIT_ASSERT_EQUAL((*diagonal)*unitVec[2],3*unitVec[2]);
[5630bd]331 CPPUNIT_ASSERT_EQUAL((*diagonal)*zeroVec,zeroVec);
332
[407782]333 CPPUNIT_ASSERT_EQUAL((*perm1)*unitVec[0],unitVec[0]);
334 CPPUNIT_ASSERT_EQUAL((*perm1)*unitVec[1],unitVec[2]);
335 CPPUNIT_ASSERT_EQUAL((*perm1)*unitVec[2],unitVec[1]);
[5630bd]336 CPPUNIT_ASSERT_EQUAL((*perm1)*zeroVec,zeroVec);
337
[407782]338 CPPUNIT_ASSERT_EQUAL((*perm2)*unitVec[0],unitVec[1]);
339 CPPUNIT_ASSERT_EQUAL((*perm2)*unitVec[1],unitVec[0]);
340 CPPUNIT_ASSERT_EQUAL((*perm2)*unitVec[2],unitVec[2]);
[5630bd]341 CPPUNIT_ASSERT_EQUAL((*perm2)*zeroVec,zeroVec);
342
[407782]343 CPPUNIT_ASSERT_EQUAL((*perm3)*unitVec[0],unitVec[1]);
344 CPPUNIT_ASSERT_EQUAL((*perm3)*unitVec[1],unitVec[2]);
345 CPPUNIT_ASSERT_EQUAL((*perm3)*unitVec[2],unitVec[0]);
[5630bd]346 CPPUNIT_ASSERT_EQUAL((*perm3)*zeroVec,zeroVec);
347
[407782]348 CPPUNIT_ASSERT_EQUAL((*perm4)*unitVec[0],unitVec[2]);
349 CPPUNIT_ASSERT_EQUAL((*perm4)*unitVec[1],unitVec[1]);
350 CPPUNIT_ASSERT_EQUAL((*perm4)*unitVec[2],unitVec[0]);
[5630bd]351 CPPUNIT_ASSERT_EQUAL((*perm4)*zeroVec,zeroVec);
352
[407782]353 CPPUNIT_ASSERT_EQUAL((*perm5)*unitVec[0],unitVec[2]);
354 CPPUNIT_ASSERT_EQUAL((*perm5)*unitVec[1],unitVec[0]);
355 CPPUNIT_ASSERT_EQUAL((*perm5)*unitVec[2],unitVec[1]);
[5630bd]356 CPPUNIT_ASSERT_EQUAL((*perm5)*zeroVec,zeroVec);
357
358 Vector t = Vector(1.,2.,3.);
359 CPPUNIT_ASSERT_EQUAL((*perm1)*t,Vector(1,3,2));
360 CPPUNIT_ASSERT_EQUAL((*perm2)*t,Vector(2,1,3));
361 CPPUNIT_ASSERT_EQUAL((*perm3)*t,Vector(3,1,2));
362 CPPUNIT_ASSERT_EQUAL((*perm4)*t,Vector(3,2,1));
363 CPPUNIT_ASSERT_EQUAL((*perm5)*t,Vector(2,3,1));
364}
Note: See TracBrowser for help on using the repository browser.