acoustic_mics_coefs()

int acoustic_mics_coefs

(

)

Evaluates B and C coefficients (please refer to the user manual for details on the theory) The actual evaluation is done in parallel on all the nodes of the domain (SPMD). The matrices are partitioned according to scalapack 2d block cyclic distribution (locRows x locCols) In each loop over the microphones, krow lines of the coefficient matrices can be loaded into memory to minimize memory usage and then written on two different files:

mbcoefs.bin.<nodeid>  B coefficients
mccoefs.bin.<nodeid>  C coefficients

where <nodeid> is the rank of the node (to avoid name-clashing on multiple process running on the same node) The same files are then read in the Matrices assembly procedure (matrices.c), basically we dictate here the overall partitioning method for the run.

See also

ac_coef_body.c

matrices.c

integrals.c

                         {


   int npsib,nphib;
   int icntr,ielmb;
   int ielmu,is;
   int mpi_nelmb;
   struct panel4* mpi_elements;
   int mpi_ncntr,mpi_nelmu;
   struct vector* mpi_cntr;

   int kode;
   int bstride;
   int cstride;

   DOUBLE theta;
   DOUBLE source,doublet,ratelet;

   int gRows,gCols;     // global dim of coefficient matrix
   int locRows,locCols; // local dim of coefficient matrix

   int izero=0;

#ifdef _ACOUSTO_DEBUG
   FILE* f;
   char fname[15];
#endif

   FILE* fbcoef = NULL;
   FILE* fccoef = NULL;
   int nrow,nrowi,iloop,nloop,irow,nio,nrlast;
   char fbname[20];
   char fcname[20];
   DOUBLE* ctmp;
   int i,j;
   int ic;
   long long int perc,io,icr,ioldcurrent,icurrent,ntotal;
   int bsymmoffset,csymmoffset;

   logger(LOG_INFO,MSG_MICSCOEF_START);

   // assigning number of psi and phi points (coincident with number of elements)
   npsib = modgeominfo->nelmb;
   nphib = modgeominfo->nelmb;


   /* copying variables into local storage (just for ease of use) */
   mpi_nelmb = modgeominfo->nelmb; 
   mpi_elements = geometry->elements;

   mpi_ncntr = modgeominfo->nmics;  // control points here are = number of microphones
   mpi_nelmu = modgeominfo->ncntr;  // control points here are = number of microphones
   mpi_cntr  = geometry->mics;


#ifdef _ACOUSTO_DEBUG
  sprintf(fname,"mics.%d%d",runinfo->myrow,runinfo->mycol);
  logger(LOG_DEBUG,"filename [%s]\n",fname);
  f = fopen(fname,"w");
  for(i=0;i<mpi_ncntr;i++){
    fprintf(f,VECTOR_PRINTF,mpi_cntr[i].x, mpi_cntr[i].y, mpi_cntr[i].z);
  }
  fclose(f);
#endif


   gRows = mpi_ncntr;
   //gCols = mpi_nelmb;
   gCols = modgeominfo->ncntr;

   micscoefs = calloc(1,sizeof(struct micscoefs_struct)); // allocating main mics coefficients struct

   // Building array descriptor
   locRows = numroc_(&gRows,&runinfo->mics_block_size,&runinfo->myrow,&izero,&runinfo->nprows);
   locCols = numroc_(&gCols,&runinfo->col_block_size,&runinfo->mycol,&izero,&runinfo->npcols);

   // nrow lines each loop
   nrow  = runinfo->krow;
   if(runinfo->krow <= 0 || runinfo->krow >= locRows) nrow = locRows;

   // ALLOCATING COEFS
   //   micscoefs->B = calloc(2*nrow*locCols,sizeof(DOUBLE)); // 2 B coefficients, nrow*partition size
   //   micscoefs->C = calloc(3*nrow*locCols,sizeof(DOUBLE)); // 3 B coefficients, nrow*partition size
   //   logger(LOG_DEBUG,"Total allocated memory for microphones coefficients=%d bytes\n",
   //                     (2*nrow*locCols*sizeof(DOUBLE) + 3*nrow*locCols*sizeof(DOUBLE)));
      micscoefs->B = calloc(2*nrow*locCols*runinfo->nsymm,sizeof(DOUBLE)); // 2 B coefficients, nrow*partition size*symmetry factor
      micscoefs->C = calloc(3*nrow*locCols*runinfo->nsymm,sizeof(DOUBLE)); // 3 B coefficients, nrow*partition size*symmetry factor
      logger(LOG_DEBUG,"Total allocated memory for microphones coefficients=%d bytes\n",
                        (2*nrow*locCols*runinfo->nsymm*sizeof(DOUBLE) + 3*nrow*locCols*runinfo->nsymm*sizeof(DOUBLE)));

      ctmp = calloc(locRows,sizeof(DOUBLE));


  if(runinfo->krow > 0){ // only if krow is used
   // setting names of files (different for each node in case of multi-process run)
   sprintf(fbname,"mbcoefs.bin.%d",rank); 
   sprintf(fcname,"mccoefs.bin.%d",rank);
   logger(LOG_DEBUG,"mics coefficients will be in files \"%s\" and \"%s\"\n", fbname,fcname);
   // opening files in binary-write mode
   fbcoef = fopen(fbname,"wb"); 
   fccoef = fopen(fcname,"wb");
  }

  // evaluating local block of coefficient
  j12(locRows,nrow,&nrlast,&nloop); // calculating start row and number of loops

  ntotal = (nloop-1) * nrow * locCols * runinfo->nsymm;
  ioldcurrent = 0;
  icurrent = 0;
  ioldcurrent = -1;

  for(iloop=0;iloop<nloop;iloop++){  // LOOP ON NLOOP
   ic = 0;
   nrowi = nrow;
   if(iloop == nloop-1) nrowi = nrlast; // setting number of rows=number of last rows on final loop
   //calculating stride for arrays
   bstride = nrowi*locCols;
   cstride = nrowi*locCols;

   for(irow=0;irow<nrowi;irow++){    // LOOP ON NROW
    i= j21(irow,iloop,nrow); // calculating control point index
    sc_l2g(i,runinfo->myrow,gRows,runinfo->nprows,runinfo->mics_block_size,&icntr);
    for(j=0;j<locCols;j++){ // LOOP ON LOCCOLS
      // getting local indices
      sc_l2g(j,runinfo->mycol,gCols,runinfo->npcols,runinfo->col_block_size,&ielmu);

      for(is=0;is<runinfo->nsymm;is++){

                  bsymmoffset = is * 2*nrowi*locCols;
                  csymmoffset = is * 3*nrowi*locCols;
                  ielmb = ielmu + is*mpi_nelmu;

                  kode = 0; // singular panels flag not applicable here

                  delaya(&theta,ielmb,&mpi_elements[ielmb],mpi_cntr[icntr]);                               // evaluating acoustic delay
                  intgba(1, 0, kode, &mpi_cntr[icntr], &mpi_elements[ielmb], &source, &doublet, &ratelet); // integrating

                  ic = irow*locCols +j;

                  /*
                  micscoefs->B[(irow*locCols+j)            ] = source;
                  micscoefs->B[(irow*locCols+j) + bstride  ] = theta/runinfo->vsound;

                  micscoefs->C[(irow*locCols+j)            ] = doublet;
                  micscoefs->C[(irow*locCols+j) + cstride  ] = ratelet/runinfo->vsound ;
                  micscoefs->C[(irow*locCols+j) + cstride*2] = theta/runinfo->vsound;
                  */
                  micscoefs->B[ic            + bsymmoffset] = source;
                  micscoefs->B[ic + bstride  + bsymmoffset] = theta/runinfo->vsound;

                  micscoefs->C[ic            + csymmoffset] = doublet;
                  micscoefs->C[ic + cstride  + csymmoffset] = ratelet/runinfo->vsound ;
                  micscoefs->C[ic + cstride*2+ csymmoffset] = theta/runinfo->vsound;
        //      ic++;

                   icurrent++;
                   io = ioldcurrent*10/ntotal;
                   icr = icurrent*10/ntotal;
                   if(icr != io){
                                 perc = icurrent*100/ntotal;
                                 logger(LOG_INFO,"Completed %lld \%% \n",perc);
                   }
                   ioldcurrent = icurrent;
      }

    } // END OF LOOP ON LOCCOLS
   } // END OF LOOP ON NROW
   if(runinfo->krow > 0){
     // append on files
     nio = fwrite(micscoefs->B,sizeof(DOUBLE),2*nrowi*locCols*runinfo->nsymm,fbcoef);
     nio = fwrite(micscoefs->C,sizeof(DOUBLE),3*nrowi*locCols*runinfo->nsymm,fccoef);
   }
  } // END OF LOOP ON NLOOP

  if(runinfo->krow > 0){
   fclose(fbcoef);
   fclose(fccoef);
  }


#ifdef _ACOUSTO_DEBUG
 if(runinfo->krow > 0){ // only if krow is used
   sprintf(fbname,"mbcoefs.bin.%d",rank); sprintf(fcname,"mccoefs.bin.%d",rank);
   fbcoef = fopen(fbname,"rb"); fccoef = fopen(fcname,"rb");
 }

 sprintf(fname,"mcoef.%d%d",runinfo->myrow,runinfo->mycol);
 f = fopen(fname,"w");
  j12(locRows,nrow,&nrlast,&nloop);
  for(iloop=0;iloop<nloop;iloop++){  // LOOP ON NLOOP
   nrowi = nrow;
   if(iloop == nloop-1) nrowi = nrlast;
   bstride = nrowi*locCols;
   cstride = nrowi*locCols;
   if(runinfo->krow > 0){ // only if krow is used
     nio = fread(micscoefs->B,sizeof(DOUBLE),2*nrowi*locCols*runinfo->nsymm,fbcoef);
     nio = fread(micscoefs->C,sizeof(DOUBLE),3*nrowi*locCols*runinfo->nsymm,fccoef);
   }
   for(irow=0;irow<nrowi;irow++){    // LOOP ON NROW
   i= j21(irow,iloop,nrow);
   sc_l2g(i,runinfo->myrow,gRows,runinfo->nprows,runinfo->mics_block_size,&icntr);
   for(j=0;j<locCols;j++){
     sc_l2g(j,runinfo->mycol,gCols,runinfo->npcols,runinfo->col_block_size,&ielmu);
     for(is=0;is<runinfo->nsymm;is++){

      ielmb = ielmu + is*mpi_nelmu;
      bsymmoffset = is * 2*nrowi*locCols;
      csymmoffset = is * 3*nrow*locCols;
      fprintf(f,"%d %d %f %f %f %f %f\n",icntr,ielmb,
              micscoefs->B[(irow*locCols+j)            + bsymmoffset],
              micscoefs->B[(irow*locCols+j) + bstride  + bsymmoffset],
              micscoefs->C[(irow*locCols+j)            + csymmoffset],
              micscoefs->C[(irow*locCols+j) + cstride  + csymmoffset],
              micscoefs->C[(irow*locCols+j) + cstride*2+ csymmoffset]);
     }
    }
   }
  }
  if(runinfo->krow > 0){ // only if krow is used
    fclose(fbcoef);
    fclose(fccoef);
  }

 fclose(f);
#endif


 if(runinfo->krow > 0){
  free(micscoefs->B);
  free(micscoefs->C);
 }

 logger(LOG_INFO,MSG_MICSCOEF_END);
 return 1;

}