Jdu spát

Former-commit-id: fa35c9818ffc5d5a57dfe1881a2a8489039dadb6

qpms/__init__.py

@@ -7,7 +7,7 @@ from sys import platform as __platform
 import warnings as __warnings
 
 try:
-    from .qpms_c import PointGroup, FinitePointGroup, FinitePointGroupElement, Particle, scatsystem_set_nthreads, ScatteringSystem, ScatteringMatrix, pitau, set_gsl_pythonic_error_handling, pgsl_ignore_error, gamma_inc, lll_reduce
+    from .qpms_c import PointGroup, FinitePointGroup, FinitePointGroupElement, Particle, scatsystem_set_nthreads, ScatteringSystem, ScatteringMatrix, pitau, set_gsl_pythonic_error_handling, pgsl_ignore_error, gamma_inc, lll_reduce, ScatteringSystemCachingMode
 except ImportError as ex:
     if __platform == "linux" or __platform == "linux2":
         if 'LD_LIBRARY_PATH' not in __os.environ.keys():

qpms/qpms_c.pyx

@@ -16,7 +16,7 @@ from .cytmatrices cimport CTMatrix, TMatrixFunction, TMatrixGenerator
 from .cymaterials cimport EpsMuGenerator
 from libc.stdlib cimport malloc, free, calloc
 import warnings
-
+import enum
 
 # Set custom GSL error handler. N.B. this is obviously not thread-safe.
 cdef char *pgsl_err_reason
@@ -331,6 +331,10 @@ cpdef void scatsystem_set_nthreads(long n):
     qpms_scatsystem_set_nthreads(n)
     return
 
+class ScatteringSystemCachingMode(enum.IntEnum):
+    NEVER = QPMS_SS_CACHE_NEVER
+    AUTO = QPMS_SS_CACHE_AUTO
+    ALWAYS = QPMS_SS_CACHE_ALWAYS
 
 cdef class ScatteringSystem:
     '''
@@ -357,7 +361,9 @@ cdef class ScatteringSystem:
         #TODO is there a way to disable the constructor outside this module?
 
     @staticmethod # We don't have any "standard" constructor for this right now
-    def create(particles, medium, FinitePointGroup sym, cdouble omega): # TODO tolerances
+    def create(particles, medium, FinitePointGroup sym, cdouble omega,
+            caching_mode = QPMS_SS_CACHE_DEFAULT): # TODO tolerances
+        caching_mode = ScatteringSystemCachingMode(caching_mode)
         # These we are going to construct
         cdef ScatteringSystem self
         cdef _ScatteringSystemAtOmega pyssw
@@ -416,6 +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)
             ss = ssw[0].ss
         finally:
             free(orig.tmg)

qpms/qpms_cdefs.pxd

@@ -550,10 +550,16 @@ cdef extern from "scatsystem.h":
         cdouble omega
         qpms_epsmu_t medium
         cdouble wavenumber
+    ctypedef enum qpms_ss_caching_mode_t:
+        QPMS_SS_CACHE_NEVER
+        QPMS_SS_CACHE_ALWAYS
+        QPMS_SS_CACHE_DEFAULT
+        QPMS_SS_CACHE_AUTO
     qpms_scatsys_at_omega_t *qpms_scatsys_apply_symmetry(const qpms_scatsys_t *orig, const qpms_finite_group_t *sym,
             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)
     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_translation_booster.c

@@ -128,6 +128,30 @@ booster_t *qpms_scatsys_translation_booster_create(
   return b;
 }
 
+int qpms_ss_create_translation_cache(qpms_scatsys_t *ss, qpms_ss_caching_mode_t m) {
+QPMS_ASSERT(ss);
+  if (ss->tbooster) {
+    QPMS_WARN("Translation cache already created?");
+    return 0;
+  }
+  switch(m) {
+    case QPMS_SS_CACHE_NEVER:
+      return 0;
+    case QPMS_SS_CACHE_AUTO:
+      QPMS_WARN("Translation operator cache heuristics not implemented, creating the cache");
+    case QPMS_SS_CACHE_ALWAYS:
+      ss->tbooster = qpms_scatsys_translation_booster_create(ss);
+      if (ss->tbooster) return 0;
+      else {
+          QPMS_WARN("Failed to create tranlation operator cache");
+          return -1;
+      }
+    default:
+      QPMS_WTF;
+  }
+  QPMS_WTF;
+}
+
 static qpms_errno_t qpms_scatsys_translation_booster_eval_bessels(
     const booster_t *b, complex double *target, complex double k // includes ref. ind. effect
     ) {

qpms/scatsystem.c

@@ -473,6 +473,8 @@ void qpms_scatsys_free(qpms_scatsys_t *ss) {
     free(ss->orbit_types);
     free(ss->saecv_sizes);
     free(ss->p_by_orbit);
+    if(ss->tbooster)
+      qpms_scatsys_translation_booster_free(ss->tbooster);
     qpms_trans_calculator_free(ss->c);
   }
   free(ss);

qpms/scatsystem.h

@@ -220,6 +220,7 @@ typedef struct qpms_scatsys_at_omega_t {
 	struct qpms_scatsysw_translation_booster *translation_cache; ///< (private) cache to speedup tranlations
 } qpms_scatsys_at_omega_t;
 
+
 /// Creates a new scatsys by applying a symmetry group onto a "proto-scatsys", copying particles if needed.
 /** In fact, it copies everything except the vswf set specs and qpms_abstract_tmatrix_t instances,
  * so keep them alive until scatsys is destroyed.
@@ -269,6 +270,25 @@ void qpms_scatsys_at_omega_free(qpms_scatsys_at_omega_t *ssw);
 qpms_scatsys_at_omega_t *qpms_scatsys_at_omega(const qpms_scatsys_t *ss,
 		complex double omega);
 
+/// Determines behaviour of qpms_ss_create_translation_cache().
+typedef enum qpms_ss_caching_mode {
+	/// Don't create the translation operator cache.
+	/** Use this if the particle positions are random. */
+	QPMS_SS_CACHE_NEVER, 
+	/// Always create the translation operator cache.
+	/** Use this for highly regular large finite structures. */
+	QPMS_SS_CACHE_ALWAYS, 
+	/// Evaluate the need of caching automatically.
+	QPMS_SS_CACHE_AUTO, 
+	QPMS_SS_CACHE_DEFAULT = QPMS_SS_CACHE_AUTO
+} qpms_ss_caching_mode_t;
+
+/// Creates some data structures for speeding up translation operator calculations.
+/** This is mostly useful for "finite lattices", where many pairs of nanoparticles
+ * share the same relative positions.
+ */
+int qpms_ss_create_translation_cache(qpms_scatsys_t *ss, qpms_ss_caching_mode_t m);
+
 /// 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.
  * 
@@ -548,7 +568,9 @@ 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);