<div dir="ltr">Dear all,<br><br>after a very useful email exchange with Jed Brown quite a while ago, I was able to find the lowest eigenvalue of a large matrix which is constructed as a tensor product. Admittedly the solution is a bit hacked, but is based on a Matrix shell and Armadillo and therefore reasonably fast. The problem seems to work well for smaller systems, but once the vectors reach a certain size, I get "out of memory" errors. I have tested the initialization of a vector of that size and multiplication by the matrix. This works fine and takes roughly 20GB of memory. There are 256 GB available, so I see no reason why the esp solvers should complain. Does anyone have an idea what goes wrong here? The error message is not very helpful and claims that a memory is requested that is way beyond any reasonable number: <i>Memory requested 18446744056529684480.</i><br>
<br>Thanks and all the best, Mathis Friesdorf<br><br><br><b>Output of the Program:</b><br><i>mathis@n180:~/localisation$ ./local_plus 27</i><br><i>System Size: 27<br>--------------------------------------------------------------------------<br>
[[30558,1],0]: A high-performance Open MPI point-to-point messaging module<br>was unable to find any relevant network interfaces:<br><br>Module: OpenFabrics (openib)<br>  Host: n180<br><br>Another transport will be used instead, although this may result in<br>
lower performance.<br>--------------------------------------------------------------------------<br>[0]PETSC ERROR: --------------------- Error Message ------------------------------------<br>[0]PETSC ERROR: Out of memory. This could be due to allocating<br>
[0]PETSC ERROR: too large an object or bleeding by not properly<br>[0]PETSC ERROR: destroying unneeded objects.<br>[0]PETSC ERROR: Memory allocated 3221286704 Memory used by process 3229827072<br>[0]PETSC ERROR: Try running with -malloc_dump or -malloc_log for info.<br>
[0]PETSC ERROR: Memory requested 18446744056529684480!<br>[0]PETSC ERROR: ------------------------------------------------------------------------<br>[0]PETSC ERROR: Petsc Release Version 3.4.4, Mar, 13, 2014 <br>[0]PETSC ERROR: See docs/changes/index.html for recent updates.<br>
[0]PETSC ERROR: See docs/faq.html for hints about trouble shooting.<br>[0]PETSC ERROR: See docs/index.html for manual pages.<br>[0]PETSC ERROR: ------------------------------------------------------------------------<br>[0]PETSC ERROR: ./local_plus on a arch-linux2-cxx-debug named n180 by mathis Wed Jun 25 12:23:01 2014<br>
[0]PETSC ERROR: Libraries linked from /home/mathis/bin_nodebug/lib<br>[0]PETSC ERROR: Configure run at Wed Jun 25 00:03:34 2014<br>[0]PETSC ERROR: Configure options PETSC_DIR=/home/mathis/petsc-3.4.4 --with-debugging=1 COPTFLAGS="-O3 -march=p4 -mtune=p4" --with-fortran=0 -with-mpi=1 --with-mpi-dir=/usr/lib/openmpi --with-clanguage=cxx --prefix=/home/mathis/bin_nodebug<br>
[0]PETSC ERROR: ------------------------------------------------------------------------<br>[0]PETSC ERROR: PetscMallocAlign() line 46 in /home/mathis/petsc-3.4.4/src/sys/memory/mal.c<br>[0]PETSC ERROR: PetscTrMallocDefault() line 189 in /home/mathis/petsc-3.4.4/src/sys/memory/mtr.c<br>
[0]PETSC ERROR: VecDuplicateVecs_Contiguous() line 62 in src/vec/contiguous.c<br>[0]PETSC ERROR: VecDuplicateVecs() line 589 in /home/mathis/petsc-3.4.4/src/vec/vec/interface/vector.c<br>[0]PETSC ERROR: EPSAllocateSolution() line 51 in src/eps/interface/mem.c<br>
[0]PETSC ERROR: EPSSetUp_KrylovSchur() line 141 in src/eps/impls/krylov/krylovschur/krylovschur.c<br>[0]PETSC ERROR: EPSSetUp() line 147 in src/eps/interface/setup.c<br>[0]PETSC ERROR: EPSSolve() line 90 in src/eps/interface/solve.c<br>
[0]PETSC ERROR: main() line 48 in local_plus.cpp<br>--------------------------------------------------------------------------<br>MPI_ABORT was invoked on rank 0 in communicator MPI_COMM_WORLD <br>with errorcode 55.<br><br>
NOTE: invoking MPI_ABORT causes Open MPI to kill all MPI processes.<br>You may or may not see output from other processes, depending on<br>exactly when Open MPI kills them.<br>--------------------------------------------------------------------------<br>
<br></i><br><b>Output of make:<br></b><i>mathis@n180:~/localisation$ make local_plus<br>mpicxx -o local_plus.o -c -Wall -Wwrite-strings -Wno-strict-aliasing -Wno-unknown-pragmas -g   -fPIC   -I/home/mathis/armadillo-4.300.8/include -lblas -llapack -L/home/mathis/armadillo-4.300.8 -O3 -larmadillo -fomit-frame-pointer -I/home/mathis/bin_nodebug/include -I/home/mathis/bin_nodebug/include -I/usr/lib/openmpi/include -I/usr/lib/openmpi/include/openmpi    -D__INSDIR__= -I/home/mathis/bin_nodebug -I/home/mathis/bin_nodebug//include -I/home/mathis/bin_nodebug/include local_plus.cpp<br>
local_plus.cpp:22:0: warning: "__FUNCT__" redefined [enabled by default]<br>In file included from /home/mathis/bin_nodebug/include/petscvec.h:10:0,<br>                 from local_plus.cpp:10:<br>/home/mathis/bin_nodebug/include/petscviewer.h:386:0: note: this is the location of the previous definition<br>
g++ -o local_plus local_plus.o -Wl,-rpath,/home/mathis/bin_nodebug//lib -L/home/mathis/bin_nodebug//lib -lslepc       -Wl,-rpath,/home/mathis/bin_nodebug/lib -L/home/mathis/bin_nodebug/lib  -lpetsc -llapack -lblas -lpthread -Wl,-rpath,/usr/lib/openmpi/lib -L/usr/lib/openmpi/lib -Wl,-rpath,/usr/lib/gcc/x86_64-linux-gnu/4.7 -L/usr/lib/gcc/x86_64-linux-gnu/4.7 -Wl,-rpath,/usr/lib/x86_64-linux-gnu -L/usr/lib/x86_64-linux-gnu -Wl,-rpath,/lib/x86_64-linux-gnu -L/lib/x86_64-linux-gnu -lmpi_cxx -lstdc++ -ldl -lmpi -lopen-rte -lopen-pal -lnsl -lutil -lgcc_s -lpthread -ldl<br>
/bin/rm -f local_plus.o<br></i><br><b>Code:<br></b><i>//Author: Mathis Friesdorf<br>//<a href="mailto:MathisFriesdorf@gmail.com">MathisFriesdorf@gmail.com</a><br><br>static char help[] = "1D chain\n";<br><br>#include <iostream><br>
#include <typeinfo><br>#include <armadillo><br>#include <petscsys.h><br>#include <petscvec.h><br>#include <petscmat.h><br>#include <slepcsys.h><br>#include <slepceps.h><br>#include <math.h><br>
#include <assert.h><br><br>PetscErrorCode BathMult(Mat H, Vec x,Vec y);<br>PetscInt       L=30,d=2,dsys;<br>PetscErrorCode ierr;<br>arma::mat hint = "1.0 0 0 0.0; 0 -1.0 2.0 0; 0 2.0 -1.0 0; 0 0 0 1.0;";<br>
<br>#define __FUNCT__ "main"<br>int main(int argc, char **argv)<br>{<br>  Mat H;<br>  EPS eps;<br>  Vec xr,xi;<br>  PetscScalar kr,ki;<br>  PetscInt j, nconv;<br> <br>  L = strtol(argv[1],NULL,10);<br>  dsys = pow(d,L);<br>
  printf("%s","System Size: ");<br>  printf("%i",L);<br>  printf("%s","\n");<br>  SlepcInitialize(&argc,&argv,(char*)0,help);<br><br>  MatCreateShell(PETSC_COMM_WORLD,dsys,dsys,dsys,dsys,NULL,&H);<br>
  MatShellSetOperation(H,MATOP_MULT,(void(*)())BathMult);<br>  ierr = MatGetVecs(H,NULL,&xr); CHKERRQ(ierr);<br>  ierr = MatGetVecs(H,NULL,&xi); CHKERRQ(ierr);<br><br>  ierr = EPSCreate(PETSC_COMM_WORLD, &eps); CHKERRQ(ierr);<br>
  ierr = EPSSetOperators(eps, H, NULL); CHKERRQ(ierr);<br>  ierr = EPSSetProblemType(eps, EPS_HEP); CHKERRQ(ierr);<br>  ierr = EPSSetWhichEigenpairs(eps,EPS_SMALLEST_REAL); CHKERRQ(ierr);<br>  ierr = EPSSetFromOptions( eps ); CHKERRQ(ierr);<br>
  ierr = EPSSolve(eps); CHKERRQ(ierr);<br>  ierr = EPSGetConverged(eps, &nconv); CHKERRQ(ierr);<br>  for (j=0; j<1; j++) {<br>    EPSGetEigenpair(eps, j, &kr, &ki, xr, xi);<br>    printf("%s","Lowest Eigenvalue: ");<br>
    PetscPrintf(PETSC_COMM_WORLD,"%9F",kr);<br>    PetscPrintf(PETSC_COMM_WORLD,"\n");<br>  }<br>  EPSDestroy(&eps);<br><br>  ierr = SlepcFinalize();<br>  return 0;<br>}<br>#undef __FUNCT__<br><br>
#define __FUNCT__ "BathMult"<br>PetscErrorCode BathMult(Mat H, Vec x, Vec y)<br>{<br>  PetscInt l;<br>  uint slice;<br>  PetscScalar *arrayin,*arrayout;<br><br>  VecGetArray(x,&arrayin);<br>  VecGetArray(y,&arrayout);<br>
  arma::cube A = arma::cube(arrayin,1,1,pow(d,L),<br>      /*copy_aux_mem*/false,/*strict*/true);<br>  arma::mat result = arma::mat(arrayout,pow(d,L),1,<br>      /*copy_aux_mem*/false,/*strict*/true);<br>  for (l=0;l<L-1;l++){<br>
    A.reshape(pow(d,L-2-l),pow(d,2),pow(d,l));<br>    result.reshape(pow(d,L-l),pow(d,l));<br>    for (slice=0;slice<A.n_slices;slice++){<br>      result.col(slice) += vectorise(A.slice(slice)*hint);<br>    }<br>  }<br>
  arrayin = A.memptr();<br>  ierr = VecRestoreArray(x,&arrayin); CHKERRQ(ierr);<br>  arrayout = result.memptr();<br>  ierr = VecRestoreArray(y,&arrayout); CHKERRQ(ierr);<br>  PetscFunctionReturn(0);<br>}<br>#undef __FUNCT__</i><b><br>
</b></div>