Commit fb3f98b5 authored by Giovanni La Mura's avatar Giovanni La Mura
Browse files

Move C6 data structure under ParticleDescriptor

parent fb3a6475
Loading
Loading
Loading
Loading
+4 −4
Original line number Diff line number Diff line
@@ -270,7 +270,7 @@ void cluster(const string& config_file, const string& data_file, const string& o
      p_output->append_line(virtual_line);
      sprintf(virtual_line, " \n");
      p_output->append_line(virtual_line);
      str(sconf, cid->c1, cid->c6);
      str(sconf, cid->c1);
      thdps(cid->c1->lm, cid->zpv);
      double exdc = sconf->exdc;
      double exri = sqrt(exdc);
@@ -797,7 +797,7 @@ int cluster_jxi488_cycle(int jxi488, ScattererConfiguration *sconf, GeometryConf
  outam0->write_to_disk(outam0_name);
  delete outam0;
#endif
  cms(cid->am, cid->c1, cid->c6);
  cms(cid->am, cid->c1);
#ifdef DEBUG_AM
  VirtualAsciiFile *outam1 = new VirtualAsciiFile();
  string outam1_name = output_path + "/c_AM1_JXI" + to_string(jxi488) + ".txt";
@@ -862,7 +862,7 @@ int cluster_jxi488_cycle(int jxi488, ScattererConfiguration *sconf, GeometryConf
#ifdef USE_NVTX
  nvtxRangePush("Average calculation");
#endif
  ztm(cid->am, cid->c1, cid->c6, cid->c9);
  ztm(cid->am, cid->c1, cid->c9);
#ifdef DEBUG_AM
  VirtualAsciiFile *outam3 = new VirtualAsciiFile();
  string outam3_name = output_path + "/c_AM3_JXI" + to_string(jxi488) + ".txt";
@@ -1056,7 +1056,7 @@ int cluster_jxi488_cycle(int jxi488, ScattererConfiguration *sconf, GeometryConf
	    }
	  }
	  // label 194
	  if (iavm == 1) crsm1(cid->vk, exri, cid->c1, cid->c6);
	  if (iavm == 1) crsm1(cid->vk, exri, cid->c1);
	  if (isam < 0) {
	    apc(cid->zpv, cid->c1->le, cid->c1->am0m, cid->c1->w, sqk, cid->gap, cid->gapp);
	    raba(cid->c1->le, cid->c1->am0m, cid->c1->w, cid->tqce, cid->tqcpe, cid->tqcs, cid->tqcps);
+2 −92
Original line number Diff line number Diff line
@@ -40,96 +40,6 @@
class ParticleDescriptor;
class mixMPI;

// /*! \brief Representation of the FORTRAN C3 blocks.
//  */
// class C3 {
// protected:
//   dcomplex *vec_tsas;
  
// public:
//   //! \brief Total forward scattering amplitude.
//   dcomplex tfsas;
//   //! \brief Total scattering amplitude.
//   dcomplex **tsas;
//   //! \brief Total geometric cross-section.
//   double gcs;
//   //! \brief Total scattering cross-section.
//   double scs;
//   //! \brief Total extinction cross-section.
//   double ecs;
//   //! \brief Total absorption cross-section.
//   double acs;

//   /*! \brief C3 instance constructor.
//    */
//   C3();

//   /*! \brief C3 instance constructor copying its contents from a preexisting object.
//    */
//   C3(const C3& rhs);

//   /*! \brief C3 instance destroyer.
//    */
//   ~C3();

// #ifdef MPI_VERSION
//   /*! \brief C3 instance constructor copying all contents off MPI broadcast from MPI process 0
//    *
//    * \param mpidata: `mixMPI *` pointer to MPI data structure.
//    */
//   C3(const mixMPI *mpidata);

//   /*! \brief send C3 instance from MPI process 0 via MPI broadcasts to all other processes
//    *
//    * \param mpidata: `mixMPI *` pointer to MPI data structure.
//    */
//   void mpibcast(const mixMPI *mpidata);
// #endif

// };

/*! \brief Representation of the FORTRAN C6 blocks.
 */
class C6 {
public:
  //! \brief LMTPO = 2 * LM + 1.
  int lmtpo;
  //! \brief QUESTION: definition?
  double *rac3j;

  /*! \brief C6 instance constructor.
   *
   * \param lmtpo: `int` QUESTION: definition?
   */
  C6(int lmtpo);

  /*! \brief C6 instance constructor copying contents from preexisting object.
   *
   * \param lmtpo: `int` QUESTION: definition?
   */
  C6(const C6& rhs);

  /*! \brief C6 instance destroyer.
   */
  ~C6();

#ifdef MPI_VERSION
  /*! \brief C6 instance constructor copying all contents off MPI broadcast from MPI process 0
   *
   * \param mpidata: `mixMPI *` pointer to MPI data structure.
   */
  C6(const mixMPI *mpidata);

  /*! \brief send C6 instance from MPI process 0 via MPI broadcasts to all other processes
   *
   * \param mpidata: `mixMPI *` pointer to MPI data structure.
   */
  void mpibcast(const mixMPI *mpidata);
#endif

};


/*! \brief Representation of the FORTRAN C9 blocks.
 */
class C9 {
@@ -223,8 +133,6 @@ class ClusterIterationData {
public:
  //! \brief Pointer to a ParticleDescriptor structure.
  ParticleDescriptor *c1;
  //! \brief Pointer to a C6 structure.
  C6 *c6;
  //! \brief Pointer to a C9 structure.
  C9 *c9;
  //! \brief Vector of geometric asymmetry factors.
@@ -565,6 +473,8 @@ public:
  double *ecscm;
  //! \brief J-vector components index matrix.
  int **ind3j;
  //! \brief J-vector boundary values. QUESTION: correct?
  double *rac3j;
  // >>> END OF SECTION NEEDED BY CLUSTER AND INCLU <<< //

  // >>> NEEDED BY INCLU <<< //
+13 −19
Original line number Diff line number Diff line
@@ -99,9 +99,8 @@ double cgev(int ipamo, int mu, int l, int m);
 *
 * \param am: `complex double **`
 * \param c1: `ParticleDescriptor *`
 * \param c6: `C6 *`
 */
void cms(dcomplex **am, ParticleDescriptor *c1, C6 *c6);
void cms(dcomplex **am, ParticleDescriptor *c1);

/*! \brief Compute orientation-averaged scattered field intensity.
 *
@@ -112,9 +111,8 @@ void cms(dcomplex **am, ParticleDescriptor *c1, C6 *c6);
 * \param vk: `double` Wave number.
 * \param exri: `double` External medium refractive index.
 * \param c1: `ParticleDescriptor *`
 * \param c6: `C6 *`
 */
void crsm1(double vk, double exri, ParticleDescriptor *c1, C6 *c6);
void crsm1(double vk, double exri, ParticleDescriptor *c1);

/*! \brief Compute the transfer vector from N2 to N1.
 *
@@ -129,11 +127,10 @@ void crsm1(double vk, double exri, ParticleDescriptor *c1, C6 *c6);
 * \param l2: `int`
 * \param m2: `int`
 * \param c1: `ParticleDescriptor *`
 * \param rac3j: `double *`
 */
dcomplex ghit_d(
	      int ihi, int ipamo, int nbl, int l1, int m1, int l2, int m2,
	      ParticleDescriptor *c1, double *rac3j
	      ParticleDescriptor *c1
);

/*! \brief Compute the transfer vector from N2 to N1.
@@ -148,12 +145,11 @@ dcomplex ghit_d(
 * \param m1: `int`
 * \param l2: `int`
 * \param m2: `int`
 * \param c1: `C1 *`
 * \param c6: `C6 *`
 * \param c1: `ParticleDescriptor *` Poiunter to a ParticleDescriptor instance.
 */
dcomplex ghit(
	      int ihi, int ipamo, int nbl, int l1, int m1, int l2, int m2,
	      ParticleDescriptor *c1, C6 *c6
	      ParticleDescriptor *c1
);

/*! \brief Compute Hankel funtion and Bessel functions.
@@ -249,9 +245,9 @@ void polar(
 *
 * \param j2: `int`
 * \param j3: `int`
 * \param c6: `C6 *` Pointer to a C6 instance.
 * \param rac3j: `double *` Vector of 3j symbols.
 */
void r3j000(int j2, int j3, C6 *c6);
void r3j000(int j2, int j3, double *rac3j);

/*! \brief Compute the 3j symbol for Clebsch-Gordan coefficients for JJ transitions.
 *
@@ -262,9 +258,9 @@ void r3j000(int j2, int j3, C6 *c6);
 * \param j3: `int`
 * \param m2: `int`
 * \param m3: `int`
 * \param c6: `C6 *`
 * \param rac3j: `double *` Vector of 3j symbols.
 */
void r3jjr(int j2, int j3, int m2, int m3, C6 *c6);
void r3jjr(int j2, int j3, int m2, int m3, double *rac3j);

/*! \brief Compute the 3j symbol for Clebsch-Gordan coefficients for JJ transitions.
 *
@@ -288,9 +284,9 @@ void r3jjr_d(int j2, int j3, int m2, int m3, double *rac3j);
 * \param j2: `int`
 * \param j3: `int`
 * \param m1: `int`
 * \param c6: `C6 *`
 * \param rac3j: `double *` Vector of 3j symbols.
 */
void r3jmr(int j1, int j2, int j3, int m1, C6 *c6);
void r3jmr(int j1, int j2, int j3, int m1, double *rac3j);

/*! \brief Compute radiation torques on a particle in Cartesian coordinates.
 *
@@ -374,9 +370,8 @@ void scr2(
 *
 * \param sconf: `ScattererConfiguration *` Pointer to scatterer configuration object.
 * \param c1: `ParticleDescriptor *` Pointer to a ParticleDescriptor instance.
 * \param c6: `C6 *` Pointer to a C6 instance.
 */
void str(ScattererConfiguration *sconf, ParticleDescriptor *c1, C6 *c6);
void str(ScattererConfiguration *sconf, ParticleDescriptor *c1);

/*! \brief Compute radiation torques on particles on a k-vector oriented system.
 *
@@ -409,9 +404,8 @@ void tqr(
 *
 * \param am: `complex double **`
 * \param c1: `ParticleDescriptor *` Pointer to a ParticleDescriptor instance.
 * \param c6: `C6 *` Pointer to a C6 instance.
 * \param c9: `C9 *` Pointer to a C9 instance.
 */
void ztm(dcomplex **am, ParticleDescriptor *c1, C6 *c6, C9 * c9);
void ztm(dcomplex **am, ParticleDescriptor *c1, C9 * c9);

#endif
+4 −5
Original line number Diff line number Diff line
@@ -53,9 +53,8 @@ void exma(dcomplex **am, dcomplex **at, ParticleDescriptor *c1);
 * \param at: `dcomplex **` TBD.
 * \param enti: `double` TBD.
 * \param c1: `ParticleDescriptor *` Pointer to a ParticleDescriptor instance.
 * \param c6: `C6 *` Pointer to a C6 instance.
 */
void incms(dcomplex **am, dcomplex **at, double enti, ParticleDescriptor *c1, C6 *c6);
void incms(dcomplex **am, dcomplex **at, double enti, ParticleDescriptor *c1);

/*! \brief C++ porting of INDME.
 *
@@ -74,13 +73,13 @@ void incms(dcomplex **am, dcomplex **at, double enti, ParticleDescriptor *c1, C6
 */
void indme(
	   int i, int npnt, int npntts, double vk, dcomplex ent, double enti,
	   dcomplex entn, int &jer, int &lcalc, dcomplex &arg, ParticleDescriptor *c1	   );
	   dcomplex entn, int &jer, int &lcalc, dcomplex &arg, ParticleDescriptor *c1
);

/*! \brief C++ porting of INSTR.
 *
 * \param sconf: `ScattererConfiguration *` Pointer to a ScattererConfiguration instance.
 * \param c1: `ParticleDescriptor *` Pointer to a ParticleDescriptor instance.
 * \param c6: `C6 *` Pointer to a C6 instance.
 */
void instr(ScattererConfiguration *sconf, ParticleDescriptor *c1, C6 *c6);
void instr(ScattererConfiguration *sconf, ParticleDescriptor *c1);
#endif // INCLUDE_INCLU_SUBS_H_
+20 −89
Original line number Diff line number Diff line
@@ -36,90 +36,6 @@
#include <mpi.h>
#endif

// C3::C3() {
//   tsas = new dcomplex*[2];
//   tsas[0] = new dcomplex[2];
//   tsas[1] = new dcomplex[2];
//   tfsas = 0.0 + 0.0 * I;
//   gcs = 0.0;
//   scs = 0.0;
//   ecs = 0.0;
//   acs = 0.0;
// }

// C3::C3(const C3& rhs) {
//   tsas = new dcomplex*[2];
//   for (int ti=0; ti<2; ti++) {
//     tsas[ti] = new dcomplex[2];
//     for (int tj=0; tj<2; tj++) tsas[ti][tj] = rhs.tsas[ti][tj];
//   }
//   tfsas = rhs.tfsas;
//   gcs = rhs.gcs;
//   scs = rhs.scs;
//   ecs = rhs.ecs;
//   acs = rhs.acs;
// }

// C3::~C3() {
//   delete[] tsas[1];
//   delete[] tsas[0];
//   delete[] tsas;
// }

// #ifdef MPI_VERSION
// C3::C3(const mixMPI *mpidata) {
//   tsas = new dcomplex*[2];
//   for (int ti=0; ti<2; ti++) {
//     tsas[ti] = new dcomplex[2];
//     MPI_Bcast(tsas[ti], 2, MPI_C_DOUBLE_COMPLEX, 0, MPI_COMM_WORLD);
//   }
//   MPI_Bcast(&tfsas, 1, MPI_C_DOUBLE_COMPLEX, 0, MPI_COMM_WORLD);
//   MPI_Bcast(&gcs, 1, MPI_DOUBLE, 0, MPI_COMM_WORLD);
//   MPI_Bcast(&scs, 1, MPI_DOUBLE, 0, MPI_COMM_WORLD);
//   MPI_Bcast(&ecs, 1, MPI_DOUBLE, 0, MPI_COMM_WORLD);
//   MPI_Bcast(&acs, 1, MPI_DOUBLE, 0, MPI_COMM_WORLD);
// }

// void C3::mpibcast(const mixMPI *mpidata) {
//   for (int ti=0; ti<2; ti++) {
//     MPI_Bcast(tsas[ti], 2, MPI_C_DOUBLE_COMPLEX, 0, MPI_COMM_WORLD);
//   }
//   MPI_Bcast(&tfsas, 1, MPI_C_DOUBLE_COMPLEX, 0, MPI_COMM_WORLD);
//   MPI_Bcast(&gcs, 1, MPI_DOUBLE, 0, MPI_COMM_WORLD);
//   MPI_Bcast(&scs, 1, MPI_DOUBLE, 0, MPI_COMM_WORLD);
//   MPI_Bcast(&ecs, 1, MPI_DOUBLE, 0, MPI_COMM_WORLD);
//   MPI_Bcast(&acs, 1, MPI_DOUBLE, 0, MPI_COMM_WORLD);
// }
// #endif

C6::C6(int _lmtpo) {
  lmtpo = _lmtpo;
  rac3j = new double[lmtpo]();
}

C6::C6(const C6& rhs) {
  lmtpo = rhs.lmtpo;
  rac3j = new double[lmtpo]();
  for (int ri = 0; ri < lmtpo; ri++) rac3j[ri] = rhs.rac3j[ri];
}

C6::~C6() {
  delete[] rac3j;
}

#ifdef MPI_VERSION
C6::C6(const mixMPI *mpidata) {
  MPI_Bcast(&lmtpo, 1, MPI_INT, 0, MPI_COMM_WORLD);
  rac3j = new double[lmtpo]();
  MPI_Bcast(rac3j, lmtpo, MPI_DOUBLE, 0, MPI_COMM_WORLD);
}

void C6::mpibcast(const mixMPI *mpidata) {
  MPI_Bcast(&lmtpo, 1, MPI_INT, 0, MPI_COMM_WORLD);
  MPI_Bcast(rac3j, lmtpo, MPI_DOUBLE, 0, MPI_COMM_WORLD);
}
#endif

C9::C9(int ndi, int _nlem, int ndit, int _nlemt) {
  gis_size_0 = ndi;
  sam_size_0 = ndit;
@@ -234,7 +150,6 @@ mixMPI::~mixMPI() {

ClusterIterationData::ClusterIterationData(GeometryConfiguration *gconf, ScattererConfiguration *sconf, const mixMPI *mpidata, const int device_count) {
  c1 = new ParticleDescriptorCluster(gconf, sconf);
  c6 = new C6(c1->lmtpo);
  const int ndi = c1->nsph * c1->nlim;
  const np_int ndit = 2 * ndi;
  c9 = new C9(ndi, c1->nlem, 2 * ndi, 2 * c1->nlem);
@@ -337,7 +252,6 @@ ClusterIterationData::ClusterIterationData(GeometryConfiguration *gconf, Scatter

ClusterIterationData::ClusterIterationData(const ClusterIterationData& rhs) {
  c1 = new ParticleDescriptorCluster(reinterpret_cast<ParticleDescriptorCluster &>(*(rhs.c1)));
  c6 = new C6(*(rhs.c6));
  const int ndi = c1->nsph * c1->nlim;
  const np_int ndit = 2 * ndi;
  c9 = new C9(*(rhs.c9));
@@ -484,7 +398,6 @@ ClusterIterationData::ClusterIterationData(const ClusterIterationData& rhs) {
#ifdef MPI_VERSION
ClusterIterationData::ClusterIterationData(const mixMPI *mpidata, const int device_count) {
  c1 = new ParticleDescriptorCluster(mpidata);
  c6 = new C6(mpidata);
  const int ndi = c1->nsph * c1->nlim;
  const np_int ndit = 2 * ndi;
  c9 = new C9(mpidata);
@@ -619,7 +532,6 @@ ClusterIterationData::ClusterIterationData(const mixMPI *mpidata, const int devi

void ClusterIterationData::mpibcast(const mixMPI *mpidata) {
  c1->mpibcast(mpidata);
  c6->mpibcast(mpidata);
  const int ndi = c1->nsph * c1->nlim;
  const np_int ndit = 2 * ndi;
  c9->mpibcast(mpidata);
@@ -703,7 +615,6 @@ ClusterIterationData::~ClusterIterationData() {
  }
  delete[] zpv;
  delete c1;
  delete c6;
  delete c9;
  delete[] gaps;
  for (int ti = 1; ti > -1; ti--) {
@@ -882,6 +793,7 @@ ParticleDescriptor::ParticleDescriptor(GeometryConfiguration *gconf, ScattererCo
  scscm = NULL;
  ecscm = NULL;
  ind3j = NULL;
  rac3j = NULL;
  // >>> NEEDED BY INCLU <<<
  rm0 = NULL;
  re0 = NULL;
@@ -991,6 +903,12 @@ ParticleDescriptor::ParticleDescriptor(const ParticleDescriptor &rhs) {
  gcsv = NULL;
  // >>> NEEDED BY CLUSTER <<<
  vintt = NULL;
  tfsas = 0.0 + I * 0.0;
  tsas = NULL;
  gcs = 0.0;
  scs = 0.0;
  ecs = 0.0;
  acs = 0.0;
  // >>> NEEDED BY CLUSTER AND INCLU <<<
  vec_am0m = NULL;
  vec_fsac = NULL;
@@ -1017,6 +935,7 @@ ParticleDescriptor::ParticleDescriptor(const ParticleDescriptor &rhs) {
  scscm = NULL;
  ecscm = NULL;
  ind3j = NULL;
  rac3j = NULL;
  // >>> NEEDED BY INCLU <<<
  rm0 = NULL;
  re0 = NULL;
@@ -1191,6 +1110,7 @@ void ParticleDescriptor::mpibcast(const mixMPI *mpidata) {
    MPI_Bcast(scscm, 2, MPI_DOUBLE, 0, MPI_COMM_WORLD);
    MPI_Bcast(ecscm, 2, MPI_DOUBLE, 0, MPI_COMM_WORLD);
    MPI_Bcast(v3j0, _nv3j, MPI_DOUBLE, 0, MPI_COMM_WORLD);
    MPI_Bcast(rac3j, _lmtpo, MPI_DOUBLE, 0, MPI_COMM_WORLD);
  }
  // >>> NEEDED BY INCLU <<< //
  if (_class_type == INCLUSION_TYPE) {
@@ -1266,6 +1186,7 @@ ParticleDescriptor::~ParticleDescriptor() {
    delete[] scscm;
    delete[] ecscm;
    delete[] ind3j;
    delete[] rac3j;
  }
  // Cluster/sphere class members, destroyed if sub-class is CLUSTER or SPHERE
  if (_class_type == CLUSTER_TYPE || _class_type == SPHERE_TYPE) {
@@ -1371,6 +1292,7 @@ ParticleDescriptorCluster::ParticleDescriptorCluster(GeometryConfiguration *gcon
  ecscm = new double[2]();
  ind3j = new int*[_lm + 1];
  for (int ii = 0; ii <= _lm; ii++) ind3j[ii] = vec_ind3j + (_lm * ii);
  rac3j = new double[_lmtpo]();
}

ParticleDescriptorCluster::ParticleDescriptorCluster(const ParticleDescriptorCluster &rhs) : ParticleDescriptor(rhs) {
@@ -1477,6 +1399,8 @@ ParticleDescriptorCluster::ParticleDescriptorCluster(const ParticleDescriptorClu
  for (int vj = 0; vj < _nv3j; vj++) v3j0[vj] = rhs.v3j0[_nv3j];
  ind3j = new int*[_lm + 1];
  for (int ii = 0; ii <= _lm; ii++) ind3j[ii] = vec_ind3j + (_lm * ii);
  rac3j = new double[_lmtpo];
  for (int ri = 0; ri < _lmtpo; ri++) rac3j[ri] = rhs.rac3j[ri];
}

#ifdef MPI_VERSION
@@ -1578,6 +1502,8 @@ ParticleDescriptorCluster::ParticleDescriptorCluster(const mixMPI *mpidata) : Pa
  MPI_Bcast(v3j0, _nv3j, MPI_DOUBLE, 0, MPI_COMM_WORLD);
  ind3j = new int*[_lm + 1];
  for (int ii = 0; ii <= _lm; ii++) ind3j[ii] = vec_ind3j + (_lm * ii);
  rac3j = new double[_lmtpo];
  MPI_Bcast(rac3j, _lmtpo, MPI_DOUBLE, 0, MPI_COMM_WORLD);
}
#endif // MPI_VERSION
// >>> End of ParticleDescriptorCluster class implementation. <<< //
@@ -1629,6 +1555,7 @@ ParticleDescriptorInclusion::ParticleDescriptorInclusion(GeometryConfiguration *
  ecscm = new double[2]();
  ind3j = new int*[_lm + 1];
  for (int ii = 0; ii <= _lm; ii++) ind3j[ii] = vec_ind3j + (_lm * ii);
  rac3j = new double[_lmtpo]();
  // Needed by INCLU
  rm0 = new dcomplex[_le]();
  re0 = new dcomplex[_le]();
@@ -1699,6 +1626,8 @@ ParticleDescriptorInclusion::ParticleDescriptorInclusion(const ParticleDescripto
  for (int vj = 0; vj < _nv3j; vj++) v3j0[vj] = rhs.v3j0[_nv3j];
  ind3j = new int*[_lm + 1];
  for (int ii = 0; ii <= _lm; ii++) ind3j[ii] = vec_ind3j + (_lm * ii);
  rac3j = new double[_lmtpo];
  for (int ri = 0; ri < _lmtpo; ri++) rac3j[ri] = rhs.rac3j[ri];
  // >>> NEEDED BY INCLU <<< //
  rm0 = new dcomplex[_le];
  re0 = new dcomplex[_le];
@@ -1776,6 +1705,8 @@ ParticleDescriptorInclusion::ParticleDescriptorInclusion(const mixMPI *mpidata)
  MPI_Bcast(v3j0, _nv3j, MPI_DOUBLE, 0, MPI_COMM_WORLD);
  ind3j = new int*[_lm + 1];
  for (int ii = 0; ii <= _lm; ii++) ind3j[ii] = vec_ind3j + (_lm * ii);
  rac3j = new double[_lmtpo];
  MPI_Bcast(rac3j, _lmtpo, MPI_DOUBLE, 0, MPI_COMM_WORLD);
  // >>> NEEDED BY INCLU <<< //
  rm0 = new dcomplex[_le];
  MPI_Bcast(rm0, _le, MPI_C_DOUBLE_COMPLEX, 0, MPI_COMM_WORLD);
Loading