Booster constructor/destructor calls etc.

Former-commit-id: 8bf3c410498ae79a7bfb22b6644b465657b62752

qpms/qpms_c.pyx

@@ -422,7 +422,7 @@ cdef class ScatteringSystem:
                 orig.p[pi].pos = p.cval().pos
                 orig.p[pi].tmatrix_id = tmindices[tm_derived_key]
             ssw = qpms_scatsys_apply_symmetry(&orig, sym.rawpointer(), omega, &QPMS_TOLERANCE_DEFAULT)
-            qpms_ss_create_translation_cache(ssw, caching_mode)
+            qpms_ss_create_translation_cache(ss, caching_mode)
             ss = ssw[0].ss
         finally:
             free(orig.tmg)
@@ -706,6 +706,16 @@ cdef class _ScatteringSystemAtOmega:
     property nirreps: 
         def __get__(self): return self.ss_pyref.nirreps
 
+    def add_translation_cache(self):
+        '''
+        Adds translation cache if the parent ScatteringSystem
+        contains the required metadata (this is not done automatically,
+        as the cache is not useful if no translation operators are to
+        be evaluated
+        '''
+        self.check()
+        qpms_ssw_create_translation_cache(self.ssw)
+
     def modeproblem_matrix_full(self):
         self.check()
         cdef size_t flen = self.ss_pyref.s[0].fecv_size

qpms/qpms_cdefs.pxd

@@ -559,7 +559,8 @@ cdef extern from "scatsystem.h":
             cdouble omega, const qpms_tolerance_spec_t *tol)
     qpms_scatsys_at_omega_t *qpms_scatsys_at_omega(const qpms_scatsys_t *ss, cdouble omega)
     void qpms_scatsys_at_omega_free(qpms_scatsys_at_omega_t *ssw)
-    int qpms_ss_create_translation_cache(qpms_scatsys_t *ss, qpms_caching_mode_t m)
+    int qpms_ss_create_translation_cache(qpms_scatsys_t *ss, qpms_ss_caching_mode_t m)
+    int qpms_ssw_create_translation_cache(qpms_scatsys_at_omega_t *ssw)
     cdouble *qpms_scatsys_irrep_pack_matrix(cdouble *target_packed,
             const cdouble *orig_full, const qpms_scatsys_t *ss, qpms_iri_t iri)
     cdouble *qpms_scatsys_irrep_unpack_matrix(cdouble *target_full, 

qpms/scatsys_private.h

@@ -0,0 +1,21 @@
+#ifndef QPMS_SCATSYS_PRIVATE_H
+#define QPMS_SCATSYS_PRIVATE_H
+
+#include "scatsystem.h"
+
+complex double *qpms_scatsysw_build_modeproblem_matrix_full_boosted(
+		complex double *target, const qpms_scatsys_at_omega_t *ssw);
+
+/// "private" destructor, called by qpms_scatsys_free()
+void qpms_scatsys_translation_booster_free(struct qpms_scatsys_translation_booster *);
+/// "private" constructor, use qpms_ss_create_translation_cache() instead.
+struct qpms_scatsys_translation_booster *qpms_scatsys_translation_booster_create(
+                const qpms_scatsys_t *ss);
+
+
+struct qpms_scatsysw_translation_booster *
+qpms_scatsysw_translation_booster_create(const qpms_scatsys_at_omega_t *ssw);
+
+/// "private" destructor, called by qpms_scatsys_at_omega_free()
+void qpms_scatsysw_translation_booster_free(struct qpms_scatsysw_translation_booster *);
+#endif 

qpms/scatsys_translation_booster.c

@@ -1,7 +1,7 @@
 // Functionality to speedup translation matrix computations in large finite arrays
 // by caching Bessel function values etc.
 #include <cblas.h>
-#include "scatsystem.h"
+#include "scatsys_private.h"
 #include "translations_inlines.h"
 #include <assert.h>
 #include <stdlib.h>
@@ -128,6 +128,15 @@ booster_t *qpms_scatsys_translation_booster_create(
   return b;
 }
 
+void qpms_scatsys_translation_booster_free(booster_t *b) {
+  if (b) {
+    free(b->bessel_offsets_r);
+    free(b->lMax_r);
+    free(b->r);
+    free(b);
+  }
+}
+
 int qpms_ss_create_translation_cache(qpms_scatsys_t *ss, qpms_ss_caching_mode_t m) {
 QPMS_ASSERT(ss);
   if (ss->tbooster) {

qpms/scatsystem.c

@@ -6,7 +6,7 @@
 #include <lapacke.h>
 #include <cblas.h>
 #include <lapacke.h>
-#include "scatsystem.h"
+#include "scatsys_private.h"
 #include "indexing.h"
 #include "vswf.h"
 #include "groups.h"
@@ -181,6 +181,7 @@ qpms_scatsys_at_omega_t *qpms_scatsys_apply_symmetry(const qpms_scatsys_t *orig,
   ss->lenscale = lenscale;
   ss->sym = sym;
   ss->medium = orig->medium;
+  ss->tbooster = NULL;
 
   // Copy the qpms_tmatrix_fuction_t from orig
   ss->tmg_count = orig->tmg_count;
@@ -521,15 +522,26 @@ qpms_scatsys_at_omega_t *qpms_scatsys_at_omega(const qpms_scatsys_t *ss,
   for (qpms_ss_tmgi_t tmgi = 0; tmgi < ss->tmg_count; ++tmgi)
     qpms_tmatrix_free(tmatrices_preop[tmgi]);
   free(tmatrices_preop);
+  ssw->translation_cache= NULL;
   return ssw;
 }
 
+int qpms_ssw_create_translation_cache(qpms_scatsys_at_omega_t *ssw) {
+  if(!ssw->translation_cache)
+    if(ssw->ss->tbooster) 
+      ssw->translation_cache = qpms_scatsysw_translation_booster_create(ssw);
+  QPMS_ASSERT(ssw->translation_cache);
+  return 0;
+}
+
 void qpms_scatsys_at_omega_free(qpms_scatsys_at_omega_t *ssw) {
   if (ssw) {
     if(ssw->tm) 
      for(qpms_ss_tmi_t i = 0; i < ssw->ss->tm_count; ++i) 
       qpms_tmatrix_free(ssw->tm[i]); 
     free(ssw->tm);
+    if(ssw->translation_cache)
+      qpms_scatsysw_translation_booster_free(ssw->translation_cache);
   }
   free(ssw);
 }

qpms/scatsystem.h

@@ -289,6 +289,10 @@ typedef enum qpms_ss_caching_mode {
  */
 int qpms_ss_create_translation_cache(qpms_scatsys_t *ss, qpms_ss_caching_mode_t m);
 
+/// Pre-calculates the "translation cache" of a qpms_scatsys_at_omega_t.
+/** If ssw->ss->tbooster is NULL, does nothing. */
+int qpms_ssw_create_translation_cache(qpms_scatsys_at_omega_t *ssw);
+
 /// Creates a "full" transformation matrix U that takes a full vector and projects it onto an symmetry adapted basis.
 /** Mostly as a reference and a debugging tool, as multiplicating these big matrices would be inefficient.
  * 
@@ -568,11 +572,6 @@ struct beyn_result_t *qpms_scatsys_finite_find_eigenmodes(
 		double res_tol ///< (default: `0.0`) TODO DOC.
 		);
 
-/// "private" destructor, called by qpms_scatsys_free()
-void qpms_scatsys_translation_booster_free(struct qpms_scatsys_translation_booster *);
-/// "private" constructor, use qpms_ss_create_translation_cache() instead.
-struct qpms_scatsys_translation_booster *qpms_scatsys_translation_booster_create(
-		const qpms_scatsys_t *ss);
 
 #if 0
 /// Searches for scattering system's eigenmodes using Beyn's algorithm.