From ebefc07e9db33c88efb1a72c51d158fa37656654 Mon Sep 17 00:00:00 2001
From: =?UTF-8?q?Marek=20Ne=C4=8Dada?= <marek@necada.org>
Date: Mon, 15 Jul 2019 14:39:08 +0300
Subject: [PATCH] Implement plane wave expansion on scatsys; untested.

Former-commit-id: d54d13d1e1994a0b9cba0860e7095063b23ff276
---
 qpms/scatsystem.c | 30 +++++++++++++++++++++++++++++-
 qpms/scatsystem.h |  7 +++++++
 2 files changed, 36 insertions(+), 1 deletion(-)

diff --git a/qpms/scatsystem.c b/qpms/scatsystem.c
index e30f94e..6ed93f9 100644
--- a/qpms/scatsystem.c
+++ b/qpms/scatsystem.c
@@ -1461,6 +1461,33 @@ complex double *qpms_scatsys_build_modeproblem_matrix_irrep_packed_parallelR(
   return target_packed;
 }
 
+
+complex double *qpms_scatsys_incident_field_vector_full(
+    complex double *target_full, const qpms_scatsys_t *ss,
+    qpms_incfield_t f,	const void *args, bool add ) {
+  QPMS_UNTESTED;
+  if (!target_full) QPMS_CRASHING_CALLOC(target_full, ss->fecv_size,
+      sizeof(complex double));
+  for(qpms_ss_pi_t pi = 0; pi < ss->p_count; ++pi) {
+    complex double *ptarget = target_full + ss->fecv_pstarts[pi];
+    const qpms_vswf_set_spec_t *bspec = qpms_ss_bspec_pi(ss, pi);
+    const cart3_t pos = ss->p[pi].pos;
+    QPMS_ENSURE_SUCCESS(f(ptarget, bspec, pos, args, add));
+  }
+  return target_full;
+}
+
+
+#if 0
+complex double *qpms_scatsys_incident_field_vector_irrep_packed(
+    complex double *target_full, const qpms_scatsys_t *ss,
+    const qpms_iri_t iri, qpms_incfield_t f,
+    const void *args, bool add) {
+  TODO;
+}
+#endif
+
+
 ccart3_t qpms_scatsys_eval_E(const qpms_scatsys_t *ss, 
     const complex double *cvf, const cart3_t where,
     const complex double k) {
@@ -1469,7 +1496,7 @@ ccart3_t qpms_scatsys_eval_E(const qpms_scatsys_t *ss,
   ccart3_t res_kc = {0,0,0}; // kahan sum compensation
 
   for (qpms_ss_pi_t pi = 0; pi < ss->p_count; ++pi) {
-    const qpms_vswf_set_spec_t *bspec = ss->tm[ss->p[pi].tmatrix_id]->spec;
+    const qpms_vswf_set_spec_t *bspec = qpms_ss_bspec_pi(ss, pi);
     const cart3_t particle_pos = ss->p[pi].pos;
     const complex double *particle_cv = cvf + ss->fecv_pstarts[pi];
 
@@ -1491,3 +1518,4 @@ ccart3_t qpms_scatsys_eval_E_irrep(const qpms_scatsys_t *ss,
   TODO;
 }
 #endif
+
diff --git a/qpms/scatsystem.h b/qpms/scatsystem.h
index 9706aac..945c305 100644
--- a/qpms/scatsystem.h
+++ b/qpms/scatsystem.h
@@ -167,6 +167,11 @@ typedef struct qpms_scatsys_t {
 	struct qpms_trans_calculator *c;
 } qpms_scatsys_t;
 
+/// Convenience function to access pi'th particle's bspec.
+static inline const qpms_vswf_set_spec_t *qpms_ss_bspec_pi(const qpms_scatsys_t *ss, qpms_ss_pi_t pi) {
+	return ss->tm[ss->p[pi].tmatrix_id]->spec;
+}
+
 /// Creates a new scatsys by applying a symmetry group, copying particles if needed.
 /** In fact, it copies everything except the vswf set specs, so keep them alive until scatsys is destroyed.
  * The following fields must be filled:
@@ -334,6 +339,7 @@ complex double *qpms_scatsys_incident_field_vector_full(
 		bool add ///< If true, add to target_full; rewrite target_full if false.
 		);
 
+#if 0
 /// Creates a (partial) incident field vector in the symmetry-adapted basis, given a function that evaluates the field expansions at points.
 /** TODO detailed doc! */
 complex double *qpms_scatsys_incident_field_vector_irrep_packed(
@@ -346,6 +352,7 @@ complex double *qpms_scatsys_incident_field_vector_irrep_packed(
 		const void *args, ///< Pointer passed as the last argument to (*field_at_point)()
 		bool add ///< If true, add to target_full; rewrite target_full if false.
 		);
+#endif
 
 /// Evaluates scattered fields (corresponding to a given excitation vector) at a given point.
 /**