UnitTests_DMRGBlock.cpp
1 static char help[] =
2  "Test code for DMRGBlock module\n";
3 
4 #include <slepceps.h>
5 #include "DMRGBlock.hpp"
6 #include "linalg_tools.hpp"
7 #include <iostream>
8 
9 PETSC_EXTERN PetscErrorCode SetRow(const Mat& A, const PetscInt& row, const std::vector<PetscInt>& idxn);
10 PETSC_EXTERN PetscErrorCode Makedir(const std::string& dir_name);
11 PETSC_EXTERN PetscErrorCode SetSz0(const Mat& Sz);
12 PETSC_EXTERN PetscErrorCode SetSp0(const Mat& Sp);
13 PETSC_EXTERN PetscErrorCode SetSz1(const Mat& Sz);
14 PETSC_EXTERN PetscErrorCode SetSp1(const Mat& Sp);
15 
16 static char hborder[] = //"************************************************************"
17  //"************************************************************";
18  "------------------------------------------------------------"
19  "------------------------------------------------------------";
20 
21 #define PrintHeader(COMM,TEXT) \
22  ierr = PetscPrintf((COMM), "%s\n%s\n%s\n", hborder, (TEXT), hborder); CHKERRQ(ierr);
23 
24 
26 PetscErrorCode Test_InitAndCopy()
27 {
28  PetscErrorCode ierr = 0;
29 
30  Block::SpinBase blk;
31  ierr = blk.Initialize(PETSC_COMM_WORLD, 2, {1.5,0.5,-0.5,-1.5}, {2,3,2,1});CHKERRQ(ierr);
32  ierr = blk.CheckSectors(); CHKERRQ(ierr);
33 
34  {
35  /* Set the entries of Sz(0) following the correct sectors */
36  ierr = SetSz0(blk.Sz(0)); CHKERRQ(ierr);
37  ierr = MatPeek(blk.Sz(0), "blk.Sz(0)"); CHKERRQ(ierr);
38  ierr = blk.MatOpCheckOperatorBlocks(OpSz, 0); CHKERRQ(ierr);
39 
40  /* Set the entries of Sp(0) following the correct sectors */
41  ierr = SetSp0(blk.Sp(0)); CHKERRQ(ierr);
42  ierr = MatPeek(blk.Sp(0), "blk.Sp(0)"); CHKERRQ(ierr);
43  ierr = blk.MatOpCheckOperatorBlocks(OpSp, 0); CHKERRQ(ierr);
44 
45  /* Set the entries of Sz(1) following the correct sectors */
46  ierr = SetSz1(blk.Sz(1)); CHKERRQ(ierr);
47  ierr = MatPeek(blk.Sz(1), "blk.Sz(1)"); CHKERRQ(ierr);
48  ierr = blk.MatOpCheckOperatorBlocks(OpSz, 1); CHKERRQ(ierr);
49 
50  /* Set the entries of Sp(1) following the correct sectors */
51  ierr = SetSp1(blk.Sp(1)); CHKERRQ(ierr);
52  ierr = MatPeek(blk.Sp(1), "blk.Sp(1)"); CHKERRQ(ierr);
53  ierr = blk.MatOpCheckOperatorBlocks(OpSp, 1); CHKERRQ(ierr);
54  }
55 
56  /* Copy to blk2 */
57  Block::SpinBase blk2 = blk;
58 
59  ierr = MatPeek(blk2.Sz(0), "blk2.Sz(0)"); CHKERRQ(ierr);
60  ierr = MatPeek(blk2.Sp(0), "blk2.Sp(0)"); CHKERRQ(ierr);
61  ierr = MatPeek(blk2.Sz(1), "blk2.Sz(1)"); CHKERRQ(ierr);
62  ierr = MatPeek(blk2.Sp(1), "blk2.Sp(1)"); CHKERRQ(ierr);
63 
64  ierr = blk2.MatOpCheckOperatorBlocks(OpSz, 0); CHKERRQ(ierr);
65  ierr = blk2.MatOpCheckOperatorBlocks(OpSp, 0); CHKERRQ(ierr);
66  ierr = blk2.MatOpCheckOperatorBlocks(OpSz, 1); CHKERRQ(ierr);
67  ierr = blk2.MatOpCheckOperatorBlocks(OpSp, 1); CHKERRQ(ierr);
68 
69  /* Call Destroy only on blk2 which should also destroy matrices in blk */
70  ierr = blk2.Destroy(); CHKERRQ(ierr);
71 
72  return ierr;
73 }
74 
76 PetscErrorCode Test_MatOpCheckOperatorBlocks()
77 {
78  PetscErrorCode ierr = 0;
79 
80  Block::SpinBase blk;
81  ierr = blk.Initialize(PETSC_COMM_WORLD, 2, {1.5,0.5,-0.5,-1.5}, {2,3,2,1});CHKERRQ(ierr);
82  ierr = blk.CheckSectors(); CHKERRQ(ierr);
83 
84  /* Set the entries of Sz(0) following the correct sectors */
85  ierr = SetRow(blk.Sz(0), 0, {0,1}); CHKERRQ(ierr);
86  ierr = SetRow(blk.Sz(0), 1, {1}); CHKERRQ(ierr);
87  ierr = SetRow(blk.Sz(0), 2, {2,3,4}); CHKERRQ(ierr);
88  ierr = SetRow(blk.Sz(0), 3, {3}); CHKERRQ(ierr);
89  ierr = SetRow(blk.Sz(0), 4, {2,4}); CHKERRQ(ierr);
90  ierr = SetRow(blk.Sz(0), 5, {5,6}); CHKERRQ(ierr);
91  ierr = SetRow(blk.Sz(0), 7, {7}); CHKERRQ(ierr);
92  ierr = MatPeek(blk.Sz(0), "blk.Sz(0)"); CHKERRQ(ierr);
93  ierr = blk.MatOpCheckOperatorBlocks(OpSz, 0); CHKERRQ(ierr);
94 
95  /* Set the entries of Sp(0) following the correct sectors */
96  ierr = SetRow(blk.Sp(0), 0, {2,3,4}); CHKERRQ(ierr);
97  ierr = SetRow(blk.Sp(0), 1, {2,4}); CHKERRQ(ierr);
98  ierr = SetRow(blk.Sp(0), 2, {5,6}); CHKERRQ(ierr);
99  ierr = SetRow(blk.Sp(0), 3, {5}); CHKERRQ(ierr);
100  ierr = SetRow(blk.Sp(0), 4, {6}); CHKERRQ(ierr);
101  ierr = SetRow(blk.Sp(0), 6, {7}); CHKERRQ(ierr);
102  ierr = MatPeek(blk.Sp(0), "blk.Sp(0)"); CHKERRQ(ierr);
103  ierr = blk.MatOpCheckOperatorBlocks(OpSp, 0); CHKERRQ(ierr);
104 
105  /* Set the entries of Sz(1) following the INCORRECT sectors */
106  ierr = SetRow(blk.Sz(1), 0, {0,1}); CHKERRQ(ierr);
107  ierr = SetRow(blk.Sz(1), 1, {1}); CHKERRQ(ierr);
108  ierr = SetRow(blk.Sz(1), 2, {2,3,4}); CHKERRQ(ierr);
109  ierr = SetRow(blk.Sz(1), 3, {3}); CHKERRQ(ierr);
110  ierr = SetRow(blk.Sz(1), 4, {2,4}); CHKERRQ(ierr);
111  ierr = SetRow(blk.Sz(1), 5, {5,6}); CHKERRQ(ierr);
112  ierr = SetRow(blk.Sz(1), 7, {1,7}); CHKERRQ(ierr); /* Error here */
113  ierr = MatPeek(blk.Sz(1), "blk.Sz(1)"); CHKERRQ(ierr);
114  ierr = blk.MatOpCheckOperatorBlocks(OpSz, 1);
115  /* Verify that error has been caught */
116  {
117  PetscInt rstart, rend;
118  MatGetOwnershipRange(blk.Sz(1), &rstart, &rend);
119  if(rstart <= 7 && 7 < rend){
120  if(ierr!=PETSC_ERR_ARG_OUTOFRANGE){
121  SETERRQ(PETSC_COMM_SELF, 1, "Failed test");
122  } else {
123  printf("%s\n",hborder);
124  printf("%s Exception caught.\n",__FUNCTION__);
125  printf("%s\n",hborder);
126  }
127  }
128  }
129  ierr = blk.Destroy(); CHKERRQ(ierr);
130  return ierr;
131 }
132 
133 PetscErrorCode Test_SavingBlocks()
134 {
135  PetscErrorCode ierr = 0;
136 
137  ierr = Makedir("trash_block_test/"); CHKERRQ(ierr);
138 
139  Block::SpinBase blk;
140  ierr = blk.Initialize(PETSC_COMM_WORLD, 2, {1.5,0.5,-0.5,-1.5}, {2,3,2,1});CHKERRQ(ierr);
141  ierr = blk.InitializeSave("trash_block_test/"); CHKERRQ(ierr);
142  ierr = blk.CheckSectors(); CHKERRQ(ierr);
143  {
144  /* Set the entries of Sz(0) following the correct sectors */
145  ierr = SetSz0(blk.Sz(0)); CHKERRQ(ierr);
146  ierr = blk.MatOpCheckOperatorBlocks(OpSz, 0); CHKERRQ(ierr);
147 
148  /* Set the entries of Sp(0) following the correct sectors */
149  ierr = SetSp0(blk.Sp(0)); CHKERRQ(ierr);
150  ierr = blk.MatOpCheckOperatorBlocks(OpSp, 0); CHKERRQ(ierr);
151 
152  /* Set the entries of Sz(1) following the correct sectors */
153  ierr = SetSz1(blk.Sz(1)); CHKERRQ(ierr);
154  ierr = blk.MatOpCheckOperatorBlocks(OpSz, 1); CHKERRQ(ierr);
155 
156  /* Set the entries of Sp(1) following the correct sectors */
157  ierr = SetSp1(blk.Sp(1)); CHKERRQ(ierr);
158  ierr = blk.MatOpCheckOperatorBlocks(OpSp, 1); CHKERRQ(ierr);
159  }
160 
161  ierr = blk.SaveAndDestroy(); CHKERRQ(ierr);
162 
163  /* Ensure that matrices are null */
164  for(PetscInt isite = 0; isite < blk.NumSites(); ++isite){
165  if(blk.Sz(isite)!=NULL) SETERRQ1(PETSC_COMM_SELF,1,"Matrix Sz(%d) must be null.", isite);
166  if(blk.Sp(isite)!=NULL) SETERRQ1(PETSC_COMM_SELF,1,"Matrix Sp(%d) must be null.", isite);
167  }
168  if(blk.H!=NULL) SETERRQ(PETSC_COMM_SELF,1,"Matrix H must be null.");
169 
170  ierr = blk.Retrieve(); CHKERRQ(ierr);
171  /* Ensure that retrieved matrices have the same entries */
172  {
173  Block::SpinBase blk2;
174  ierr = blk2.Initialize(PETSC_COMM_WORLD, 2, {1.5,0.5,-0.5,-1.5}, {2,3,2,1});CHKERRQ(ierr);
175  ierr = blk2.InitializeSave("trash_block_test/"); CHKERRQ(ierr);
176  ierr = blk2.CheckSectors(); CHKERRQ(ierr);
177  PetscBool flg;
178  ierr = SetSz0(blk2.Sz(0)); CHKERRQ(ierr);
179  ierr = SetSp0(blk2.Sp(0)); CHKERRQ(ierr);
180  ierr = SetSz1(blk2.Sz(1)); CHKERRQ(ierr);
181  ierr = SetSp1(blk2.Sp(1)); CHKERRQ(ierr);
182  ierr = blk2.AssembleOperators(); CHKERRQ(ierr);
183 
184  ierr = MatEqual(blk2.Sz(0), blk.Sz(0), &flg); CHKERRQ(ierr);
185  if(!flg) SETERRQ(PETSC_COMM_WORLD,1,"Wrong retrieval of matrix Sz0");
186  ierr = MatEqual(blk2.Sp(0), blk.Sp(0), &flg); CHKERRQ(ierr);
187  if(!flg) SETERRQ(PETSC_COMM_WORLD,1,"Wrong retrieval of matrix Sp0");
188  ierr = MatEqual(blk2.Sz(1), blk.Sz(1), &flg); CHKERRQ(ierr);
189  if(!flg) SETERRQ(PETSC_COMM_WORLD,1,"Wrong retrieval of matrix Sz1");
190  ierr = MatEqual(blk2.Sp(1), blk.Sp(1), &flg); CHKERRQ(ierr);
191  if(!flg) SETERRQ(PETSC_COMM_WORLD,1,"Wrong retrieval of matrix Sp1");
192 
193  ierr = blk2.Destroy(); CHKERRQ(ierr);
194  }
195 
196  ierr = blk.Destroy(); CHKERRQ(ierr);
197 
198  return(0);
199 }
200 
201 
202 int main(int argc, char **argv)
203 {
204  PetscErrorCode ierr = 0;
205  PetscMPIInt nprocs, rank;
206  MPI_Comm& comm = PETSC_COMM_WORLD;
207 
208  /* Initialize MPI */
209  ierr = SlepcInitialize(&argc, &argv, (char*)0, help); CHKERRQ(ierr);
210  ierr = MPI_Comm_size(comm, &nprocs); CHKERRQ(ierr);
211  ierr = MPI_Comm_rank(comm, &rank); CHKERRQ(ierr);
212 
213  PrintHeader(comm, "Test 01: Test_InitAndCopy");
214  ierr = Test_InitAndCopy(); CHKERRQ(ierr);
215  PrintHeader(comm, "Test 02: Test_SavingBlocks");
216  ierr = Test_SavingBlocks(); CHKERRQ(ierr);
217  PrintHeader(comm, "Test 03: Test_MatOpCheckOperatorBlocks");
218  ierr = Test_MatOpCheckOperatorBlocks(); CHKERRQ(ierr);
219 
220  ierr = SlepcFinalize(); CHKERRQ(ierr);
221  return ierr;
222 }
Block::SpinBase::InitializeSave
PetscErrorCode InitializeSave(const std::string &save_dir_in)
Initializes the writing of the block matrices to file.
Definition: DMRGBlock.cpp:835
Makedir
PETSC_EXTERN PetscErrorCode Makedir(const std::string &dir_name)
Creates a directory named dir_name.
Definition: MiscTools.cpp:306
Block::SpinBase
Base class for the implementation of a block of spin sites.
Definition: DMRGBlock.hpp:79
OpSp
operator
Definition: DMRGBlock.hpp:25
Block::SpinBase::AssembleOperators
PetscErrorCode AssembleOperators()
Ensures that all operators are assembled.
Definition: DMRGBlock.cpp:825
Block::SpinBase::H
Mat H
Matrix representation of the Hamiltonian operator.
Definition: DMRGBlock.hpp:350
Block::SpinBase::CheckSectors
PetscErrorCode CheckSectors() const
Checks whether sector indexing was done properly.
Definition: DMRGBlock.cpp:429
Block::SpinBase::Retrieve
PetscErrorCode Retrieve()
Retrieve all the matrix operators that were written to file by SaveAndDestroy()
Definition: DMRGBlock.cpp:1048
Block::SpinBase::Sz
Mat Sz(const PetscInt &Isite) const
Returns the matrix pointer to the operator at site Isite.
Definition: DMRGBlock.hpp:353
Block::SpinBase::SaveAndDestroy
PetscErrorCode SaveAndDestroy()
Save all the matrix operators to file and destroy the current storage.
Definition: DMRGBlock.cpp:1014
Block::SpinBase::Destroy
PetscErrorCode Destroy()
Destroys all operator matrices and frees memory.
Definition: DMRGBlock.cpp:655
OpSz
operator
Definition: DMRGBlock.hpp:24
Block::SpinBase::MatOpCheckOperatorBlocks
PetscErrorCode MatOpCheckOperatorBlocks(const Op_t &op_t, const PetscInt &isite) const
Checks the block indexing in the matrix operator op_t on site isite.
Definition: DMRGBlock.cpp:498
Block::SpinBase::Initialize
PetscErrorCode Initialize(const MPI_Comm &comm_in)
Initializes block object&#39;s MPI attributes.
Definition: DMRGBlock.cpp:31
Block::SpinBase::NumSites
PetscInt NumSites() const
Gets the number of sites that are currently initialized.
Definition: DMRGBlock.hpp:344
Block::SpinBase::Sp
Mat Sp(const PetscInt &Isite) const
Returns the matrix pointer to the operator at site Isite.
Definition: DMRGBlock.hpp:359
This site was generated by Sphinx using Doxygen with a customized theme from doxygen-bootstrapped.