Commit 05646e46 authored by D'azevedo, Ed's avatar D'azevedo, Ed
Browse files

initial checkin

parent bddccd46
Loading
Loading
Loading
Loading
+1 −0
Original line number Diff line number Diff line
@@ -9,6 +9,7 @@ OBJS=\
	setup_matrix.o \
	setup_vbatch.o \
	unsetup_vbatch.o \
	setup_nC.o \
	setup_sparse_batch.o \
        unsetup_sparse_batch.o \
        apply_Htarget_sparse.o \

src/setup_nC.c

0 → 100644
+147 −0
Original line number Diff line number Diff line
#include <stdlib.h>
#include <assert.h>
#include <inttypes.h>
#include <math.h>

#include "dmrg_vbatch.h"
#include "setup_nC.h"


#ifndef MIN
#define MIN(x,y)  (((x) < (y))?(x):(y))
#endif


#define index3(ipatch,jpatch,ioperator) \
   (  (  ((ipatch)-1) + (((jpatch)-1)*(npatches))  ) + \
         ( ((ioperator)-1)*((npatches)*(npatches)) )   )

#define gnnz_A(ipatch,jpatch,ioperator)  gnnz_A_[index3(ipatch,jpatch,ioperator)]

#define gnnz_B(ipatch,jpatch,ioperator) gnnz_B_[ index3(ipatch,jpatch,ioperator)]

#define Amatrix(ipatch,jpatch,ioperator) Amatrix_[  index3(ipatch,jpatch,ioperator) ]
#define Bmatrix(ipatch,jpatch,ioperator) Bmatrix_[  index3(ipatch,jpatch,ioperator) ]

#define ld_Amatrix(ipatch,jpatch,ioperator) ld_Amatrix_[  index3(ipatch,jpatch,ioperator) ]
#define ld_Bmatrix(ipatch,jpatch,ioperator) ld_Bmatrix_[  index3(ipatch,jpatch,ioperator) ]



#define nC(i,j)  nC_[ indx2f(i,j,npatches)]

#define left_patch_size(ipatch) left_patch_size_[ (ipatch)-1]
#define right_patch_size(ipatch) right_patch_size_[ (ipatch)-1]


void setup_nC( int npatches,
               int noperator,

               /*
                * ---------
                * intent(in)
                * ---------
                */
               FpType **Amatrix_, int *ld_Amatrix_,
               FpType **Bmatrix_, int *ld_Bmatrix_,
               int left_patch_size_[], /* sizes of left patches (INPUT) */
               int right_patch_size_[], /* sizes of right patches (INPUT) */

               /*
                * --------------------------------
                * intent(out) but assume storage already allocated
                * --------------------------------
                */
               IntegerType nC_[],
               IntegerType gnnz_A_[],
               IntegerType gnnz_B_[]
             )

{
   int ipatch = 0;
   int jpatch = 0;
   int ioperator = 0;

   for(jpatch=1; jpatch <= npatches; jpatch++) {
   for(ipatch=1; ipatch <= npatches; ipatch++) {
     nC(ipatch,jpatch) = 0;
     };
     };



   for(ioperator=1; ioperator <= noperator; ioperator++) {
   for(jpatch=1; jpatch <= npatches; jpatch++) {
   for(ipatch=1; ipatch <= npatches; ipatch++) {

#define Amat(i,j)  Amat_[indx2f(i,j,ld_Amat)]
#define Bmat(i,j)  Bmat_[indx2f(i,j,ld_Bmat)]
      FpType *Amat_ = Amatrix(ipatch,jpatch,ioperator);
      FpType *Bmat_ = Bmatrix(ipatch,jpatch,ioperator);
      int ld_Amat = ld_Amatrix(ipatch,jpatch,ioperator);
      int ld_Bmat = ld_Bmatrix(ipatch,jpatch,ioperator);

      int nrowA = left_patch_size(ipatch);
      int ncolA = left_patch_size(jpatch);
 
      int nrowB = right_patch_size(ipatch);
      int ncolB = right_patch_size(jpatch);
      
      int nnz_Amat = 0;
      int nnz_Bmat = 0;

      /*
       --------------------------------------------------------------
       Note: treat a NULL pointer as pointer to a matrix of all zeros
       --------------------------------------------------------------
       */
      if (Amat_ == NULL) {
          nnz_Amat = 0;
          }
      else {
          int ia = 0;
          int ja = 0;
          for(ja=1; ja <= ncolA; ja++) {
          for(ia=1; ia <= nrowA; ia++) {
               nnz_Amat += ((Amat(ia,ja) == 0)? 0 : 1);
               };
               };
          };
 
      if (Bmat_ == NULL) {
         nnz_Bmat = 0;
         }
      else {
         int ib = 0;
         int jb = 0;
         for(jb=1; jb <= ncolB; jb++) {
         for(ib=1; ib <= nrowB; ib++) {
             nnz_Bmat += ((Bmat(ib,jb) == 0) ? 0 : 1);
             };
             };
          };
      gnnz_A(ipatch,jpatch,ioperator) = nnz_Amat;
      gnnz_B(ipatch,jpatch,ioperator) = nnz_Bmat;

      int is_zero_Amat = (nnz_Amat == 0);
      int is_zero_Bmat = (nnz_Bmat == 0);
      if (is_zero_Amat || is_zero_Bmat) {
         /*
          ---------------------------------------------------------
          ignore this operator since kron(Amat,Bmat) is zero matrix
          ---------------------------------------------------------
          */
         }
      else {
         /*
          ------------------------------------------
          pair of non-zero Amat(k), Bmat(k) matrices
          ------------------------------------------
          */
         nC(ipatch,jpatch) = nC(ipatch,jpatch) + 1;
         };
     }; /* end for jpatch */
     }; /* end for ipatch */
     }; /* end for ioperator */

}

src/setup_nC.h

0 → 100644
+44 −0
Original line number Diff line number Diff line
#ifndef SETUP_NC_H
#define SETUP_NC_H 1

#include "dmrg_types.h"


#ifdef __cplusplus
extern "C" {
#endif

extern
void setup_nC( int npatches,
               int noperator,
               /*
                * ---------
                * intent(in)
                * ---------
                */
               FpType **Amatrix_, int *ld_Amatrix_,
               FpType **Bmatrix_, int *ld_Bmatrix_,
               int left_patch_size_[], /* sizes of left patches (INPUT) */
               int right_patch_size_[], /* sizes of right patches (INPUT) */



               /*
                * --------------------------------
                * intent(out) but assume storage already allocated
                * --------------------------------
                */
               IntegerType nC_[],
               IntegerType gnnz_A_[],
               IntegerType gnnz_B_[]
             );


#ifdef __cplusplus
}
#endif



#endif