Prepare ParticleDescriptor base class and SPHERE specialization

};

/*! \brief Basic data structure describing the particle model and its interaction with fields.

 *

 * This class forms a base of the data structure collections that are used by the

 * the code to handle different scattering problems. The class implements real

 * members that are used by all code sections.

 */

class ParticleDescriptor {

protected:

  // >>> COMMON TO ALL DESCRIPTOR TYPES <<< //

  //! \brief Number of spheres composing the model.

  int _nsph;

  //! \brief Maximum internal field expansion order.

  int _li;

  //! \brief Number of different sphere types.

  int _num_configurations;

  //! \brief Total number of layers from all sphere types.

  int _num_layers;

  //! \brief Contiguous space for RMI.

  dcomplex *vec_rmi;

  //! \brief Contiguous space for REI.

  dcomplex *vec_rei;

  //! \brief Contiguous space for W.

  dcomplex *vec_w;

  //! \brief Contiguous space for RC.

  double *vec_rc;

  // >>> END OF SECTION COMMON TO ALL DESCRIPTOR TYPES <<< //

  // >>> NEEDED BY SPHERE AND CLUSTER <<< //

  //! \brief Contiguous space for SAS.

  dcomplex *vec_sas;

  //! \brief Contiguous space for VINTS.

  dcomplex *vec_vints;

  // >>> END OF SECTION NEEDED BY SPHERE AND CLUSTER <<< //

  // >>> NEEDED BY CLUSTER AND INCLU <<<

  //! \brief Maximum external field expansion order.

  int _le;

  //! \brief Maximum field expansion order.

  int _lm;

  //! \brief Number of rows of the Transition Matrix.

  int _nlmmt;

  //! \brief Contiguous space for AM0M

  dcomplex *vec_am0m;

  //! \brief Contiguous space for FSAC

  dcomplex *vec_fsac;

  //! \brief Contiguous space for SAC

  dcomplex *vec_sac;

  //! \brief Contiguous space for FSACM

  dcomplex *vec_fsacm;

  // >>> END OF SECTION NEEDED BY CLUSTER AND INCLU <<<

public:

  // >>> COMMON TO ALL DESCRIPTOR TYPES <<< //

  //! \brief Read-only view of number of spheres composing the model.

  const int &nsph = _nsph;

  //! \brief Read-only view of maximum internal field expansion order.

  const int &li = _li;

  //! \brief Read-only view of number of different sphere types.

  const int &num_configurations = _num_configurations;

  //! \brief Read-only view of total number of layers from all sphere types.

  const int &num_layers = _num_layers;

  //! \brief Vector of number of layers in type.

  int *num_layers_in_type;

  //! \brief Matrix of inverse radial M coefficients.

  dcomplex **rmi;

  //! \brief Matrix of inverse radial E coefficients.

  dcomplex **rei;

  //! \brief Matrix of W coefficients.

  dcomplex **w;

  //! \brief Vector of intensity components.

  dcomplex *vint;

  //! \brief Vector of sphere Cartesian X coordinates.

  double *rxx;

  //! \brief Vector of sphere Cartesian Y coordinates.

  double *ryy;

  //! \brief Vector of sphere Cartesian Z coordinates.

  double *rzz;

  //! \brief Vector of sphere radii.

  double *ros;

  //! \brief Matrix of sphere fractional transition radii.

  double **rc;

  //! \brief Vector of sphere type identifiers.

  int *iog;

  //! \brief Vector of number of layers in sphere type.

  int *nshl;

  // >>> END OF SECTION COMMON TO ALL DESCRIPTOR TYPES <<< //

  // >>> NEEDED BY SPHERE AND CLUSTER <<< //

  //! \brief Tensor of sphere scattering amplitudes.

  dcomplex ***sas;

  //! \brief Array of vectors of intensity components for single spheres.

  dcomplex **vints;

  //! \brief Vector of sphere forward scattering amplitudes.

  dcomplex *fsas;

  //! \brief Vector of scattering cross-sections for spheres.

  double *sscs;

  //! \brief Vector of extinction cross-sections for spheres.

  double *sexs;

  //! \brief Vector of absorption cross-sections for spheres.

  double *sabs;

  //! \brief Vector of scattering efficiencies for spheres.

  double *sqscs;

  //! \brief Vector of extinction efficiencies for spheres.

  double *sqexs;

  //! \brief Vector of absorption efficiencies for spheres.

  double *sqabs;

  //! \brief Vector of geometric cross-sections for spheres.

  double *gcsv;

  // >>> END OF SECTION NEEDED BY SPHERE AND CLUSTER <<<

  // >>> NEEDED BY CLUSTER <<<

  dcomplex *vintt;

  // >>> END OF SECTION NEEDED BY CLUSTER <<<

  // >>> NEEDED BY CLUSTER AND INCLU <<<

  //! \brief Read-only view of maximum external field expansion order.

  const int& le = _le;

  //! \brief Read-only view of maximum field expansion order.

  const int& lm = _lm;

  //! \brief Read-only view of number of rows of the Transition Matrix.

  const int& nlmmt = _nlmmt;

  dcomplex *vh;

  dcomplex *vj0;

  dcomplex *vyhj;

  dcomplex *vyj0;

  dcomplex vj;

  dcomplex **am0m;

  dcomplex **fsac;

  dcomplex **sac;

  dcomplex **fsacm;

  dcomplex *vintm;

  dcomplex *scscp;

  dcomplex *ecscp;

  dcomplex *scscpm;

  dcomplex *ecscpm;

  double *v3j0;

  double *scsc;

  double *ecsc;

  double *scscm;

  double *ecscm;

  int *ind3j;

  // >>> END OF SECTION NEEDED BY CLUSTER AND INCLU <<<

  /*! \brief ParticleDescriptor instance constructor.

   *

   * \param gconf: `const GeometryConfiguration *` Pointer to GeometryConfiguration instance.

   * \param sconf: `const ScattererConfiguration *` Pointer to ScattererConfiguration instance.

   */

  ParticleDescriptor(GeometryConfiguration *gconf, ScattererConfiguration *sconf);

  /*! \brief ParticleDescriptor instance destroyer.

   */

  ~ParticleDescriptor();

  /*! \brief Interface function to return the descriptor type as string.

   *

   * \return descriptor_type: `string` The descriptor type name.

   */

  virtual std::string get_descriptor_type() { return "base descriptor"; }

};

/*! \brief The data structure describing a particle model made by a cluster of spheres.

 *

 * This class is used to solve the problem of a cluster of spherical particles.

 * It replaces the C1 and C1_AddOns class implementation of the C1 FORTRAN common

 * block.

 */

class ParticleDescriptorCluster: public ParticleDescriptor {

public:

  /*! \brief ParticleDescriptorCluster instance constructor.

   *

   * \param gconf: `const GeometryConfiguration *` Pointer to GeometryConfiguration instance.

   * \param sconf: `const ScattererConfiguration *` Pointer to ScattererConfiguration instance.

   */

  ParticleDescriptorCluster(GeometryConfiguration *gconf, ScattererConfiguration *sconf);

  /*! \brief ParticleDescriptorSphere instance destroyer.

   */

  ~ParticleDescriptorCluster();

  /*! \brief Interface function to return the descriptor type as string.

   *

   * \return descriptor_type: `string` The descriptor type name.

   */

  std::string get_descriptor_type() override { return "cluster descriptor"; }

};

/*! \brief The data structure describing a spherical particle model.

 *

 * This class is used to solve the problem of a single spherical particle. It

 * replaces the old C1 class implementation of the corresponding FORTRAN common

 * block.

 */

class ParticleDescriptorSphere: public ParticleDescriptor {

public:

  /*! \brief ParticleDescriptorSphere instance constructor.

   *

   * \param gconf: `const GeometryConfiguration *` Pointer to GeometryConfiguration instance.

   * \param sconf: `const ScattererConfiguration *` Pointer to ScattererConfiguration instance.

   */

  ParticleDescriptorSphere(GeometryConfiguration *gconf, ScattererConfiguration *sconf);

  /*! \brief ParticleDescriptorSphere instance destroyer.

   */

  ~ParticleDescriptorSphere();

  /*! \brief Interface function to return the descriptor type as string.

   *

   * \return descriptor_type: `string` The descriptor type name.

   */

  std::string get_descriptor_type() override { return "sphere descriptor"; }

};

/*! \brief A data structure representing the angles to be evaluated in the problem.

 *

 */

  delete[] cmul;

}

// >>> ParticleDescriptor class implementation. <<< //

ParticleDescriptor::ParticleDescriptor(GeometryConfiguration *gconf, ScattererConfiguration *sconf) {

  _nsph = gconf->number_of_spheres;

  _li = (nsph == 1) ? gconf->l_max : gconf->li;

  _le = (nsph == 1) ? 0 : gconf->le;

  _lm = (_le > _li) ? _le : _li;

  _nlmmt = (nsph == 1) ? 2 * li * (li + 2) : 2 * _le * (_le + 2);

  _num_configurations = sconf->configurations;

  _num_layers = (sconf->use_external_sphere) ? 1 : 0;

  for (int nli = 0; nli < num_configurations; nli++) _num_layers += sconf->get_nshl(nli);

  num_layers_in_type = new int[num_configurations]();

  vec_rmi = new dcomplex[li * nsph]();

  rmi = new dcomplex*[li];

  vec_rei = new dcomplex[li * nsph]();

  rei = new dcomplex*[li];

  for (int ri = 0; ri < li; ri++) {

    rmi[ri] = vec_rmi + (nsph * ri);

    rei[ri] = vec_rei + (nsph * ri);

  }

  vec_w = new dcomplex[nllt * 4]();

  w = new dcomplex*[nllt];

  for (int wi = 0; wi < nllt; wi++) w[wi] = vec_w + (4 * wi);

  vec_rc = new double[num_layers]();

  rc = new double*[num_configurations];

  int last_layer_index = 0, cur_layer_index = 0;

  for (int rci = 0; rci < num_configurations; rci++) {

    rc[rci] = vec_rc + last_layer_index;

    last_layer_index += sconf->get_nshl(rci);

    if (rci == 0 && sconf->use_external_sphere) last_layer_index++;

    num_layers_in_type[rci] = last_layer_index - cur_layer_index;

    for (int rcj = cur_layer_index; rcj < last_layer_index; rcj++) {

      rc[rci][rcj] = sconf->get_rcf(rci, rcj);

    }

    cur_layer_index = last_layer_index;

  }

  vint = new dcomplex[16]();

  rxx = new double[nsph]();

  ryy = new double[nsph]();

  rzz = new double[nsph]();

  iog = new int[nsph]();

  for (int si = 0; si < nsph; si++) {

    rxx[si] = gconf->get_sph_x(si);

    ryy[si] = gconf->get_sph_y(si);

    rzz[si] = gconf->get_sph_z(si);

    iog[si] = sconf->get_iog(si);

  }

  ros = new double[num_configurations]();

  nshl = new int[num_configurations]();

  for (int ci = 0; ci < num_configurations; ci++) {

    ros[ci] = sconf->get_radius(ci);

    nshl[ci] = sconf->get_nshl(ci);

  }

  // >>> NEEDED BY SPHERE AND CLUSTER <<<

  sas = NULL;

  vints = NULL;

  fsas = NULL;

  sscs = NULL;

  sexs = NULL;

  sabs = NULL;

  sqscs = NULL;

  sqexs = NULL;

  sqabs = NULL;

  gcsv = NULL;

  vec_am0m = NULL;

  vec_fsac = NULL;

  vec_sac = NULL;

  vec_fsacm = NULL;

  // >>> NEEDED BY CLUSTER <<<

  vintt = NULL;

  // >>> NEEDED BY CLUSTER AND INCLU <<<

  vh = NULL;

  vj0 = NULL;

  vyhj = NULL;

  vyj0 = NULL;

  vj = 0.0 + I * 0.0;

  am0m = NULL;

  fsac = NULL;

  sac = NULL;

  fsacm = NULL;

  vintm = NULL;

  scscp = NULL;

  ecscp = NULL;

  scscpm = NULL;

  ecscpm = NULL;

  v3j0 = NULL;

  scsc = NULL;

  ecsc = NULL;

  scscm = NULL;

  ecscm = NULL;

  ind3j = NULL;

}

ParticleDescriptor::~ParticleDescriptor() {

  delete[] nshl;

  delete[] ros;

  delete[] iog;

  delete[] rzz;

  delete[] ryy;

  delete[] rxx;

  delete[] vint;

  delete[] rc;

  delete[] vec_rc;

  delete[] w;

  delete[] vec_w;

  delete[] rei;

  delete[] vec_rei;

  delete[] rmi;

  delete[] vec_rmi;

  delete[] num_layers_in_type;

}

// >>> End of ParticleDescriptor class implementation. <<< //

// >>> ParticleDescriptorSphere class implementation. <<< //

ParticleDescriptorSphere::ParticleDescriptorSphere(GeometryConfiguration *gconf, ScattererConfiguration *sconf) : ParticleDescriptor(gconf, sconf) {

  vec_sas = new dcomplex[4 * nsph]();

  vec_vints = new dcomplex[16 * nsph]();

  fsas = new dcomplex[nsph];

  sscs = new double[nsph]();

  sexs = new double[nsph]();

  sabs = new double[nsph]();

  sqscs = new double[nsph]();

  sqexs = new double[nsph]();

  sqabs = new double[nsph]();

  gcsv = new double[nsph]();

  sas = new dcomplex**[nsph];

  vints = new dcomplex*[nsph];

  for (int vi = 0; vi < nsph; vi++) {

    vints[vi] = vec_vints + (16 * nsph);

    sas[vi] = new dcomplex*[2];

    sas[vi][0] = vec_sas + (4 * vi);

    sas[vi][1] = vec_sas + (4 * vi) + 2;

  }

}

ParticleDescriptorSphere::~ParticleDescriptorSphere() {

  for (int vi = 0; vi < nsph; vi++) {

    delete[] sas[vi];

  }

  delete[] vints;

  delete[] sas;

  delete[] gcsv;

  delete[] sqabs;

  delete[] sqexs;

  delete[] sqscs;

  delete[] sabs;

  delete[] sexs;

  delete[] sscs;

  delete[] fsas;

  delete[] vec_vints;

  delete[] vec_sas;

}

// >>> End of ParticleDescriptorSphere class implementation. <<< //

// >>> ScatteringAngles class implementation. <<< //

ScatteringAngles::ScatteringAngles(GeometryConfiguration *gconf) {

  int isam = gconf->isam;

  _th = gconf->in_theta_start;

  MPI_Bcast(&_nkks, 1, MPI_INT, 0, MPI_COMM_WORLD);

}

#endif

// >>> End of ScatteringAngles class implementation. <<< //