Changeset 502c7b8 in sasmodels for sasmodels/kernel_iq.c

Timestamp:

Aug 7, 2018 5:14:59 PM (6 years ago)

Author:

Paul Kienzle <pkienzle@…>

Branches:

master, core_shell_microgels, magnetic_model, ticket-1257-vesicle-product, ticket_1156, ticket_1265_superball, ticket_822_more_unit_tests

Children:

7b0abf8

Parents:

cdd676e (diff), 01c8d9e (diff)
Note: this is a merge changeset, the changes displayed below correspond to the merge itself.
Use the (diff) links above to see all the changes relative to each parent.

Message:

Merge remote-tracking branch 'greg/beta_approx' into beta_approx

File:

• sasmodels/kernel_iq.c (modified) (8 diffs)

sasmodels/kernel_iq.c

@@ -36 +36 @@
 //  PROJECTION : equirectangular=1, sinusoidal=2
 //      see explore/jitter.py for definitions.
+
 
 #ifndef _PAR_BLOCK_ // protected block so we can include this code twice.
@@ -270 +271 @@
 #endif // _QABC_SECTION
 
-
 // ==================== KERNEL CODE ========================
-
 kernel
 void KERNEL_NAME(
@@ -339 +338 @@
   // The code differs slightly between opencl and dll since opencl is only
   // seeing one q value (stored in the variable "this_result") while the dll
-  // version must loop over all q.
+  // version must loop over all q
+  #if BETA
+  double *beta_result = &(result[nq+1]); // = result + nq + 1
+  double weight_norm = (pd_start == 0 ? 0.0 : result[2*nq+1]);
+  #endif
   #ifdef USE_OPENCL
     double pd_norm = (pd_start == 0 ? 0.0 : result[nq]);
     double this_result = (pd_start == 0 ? 0.0 : result[q_index]);
+    #if BETA
+      double this_beta_result = (pd_start == 0 ? 0.0 : result[nq + q_index];
   #else // !USE_OPENCL
-    double pd_norm = (pd_start == 0 ? 0.0 : result[nq]);
+    double pd_norm = (pd_start == 0 ? 0.0 : result[nq]); 
     if (pd_start == 0) {
       #ifdef USE_OPENMP
@@ -350 +355 @@
       #endif
       for (int q_index=0; q_index < nq; q_index++) result[q_index] = 0.0;
+      #if BETA
+      for (int q_index=0; q_index < nq; q_index++) beta_result[q_index] = 0.0;
+      #endif
     }
     //if (q_index==0) printf("start %d %g %g\n", pd_start, pd_norm, result[0]);
 #endif // !USE_OPENCL
+
+
+
 
 
@@ -442 +453 @@
 // inner loop and defines the macros that use them.
 
+
 #if defined(CALL_IQ)
   // unoriented 1D
   double qk;
-  #define FETCH_Q() do { qk = q[q_index]; } while (0)
-  #define BUILD_ROTATION() do {} while(0)
-  #define APPLY_ROTATION() do {} while(0)
-  #define CALL_KERNEL() CALL_IQ(qk, local_values.table)
+  #if BETA == 0
+    #define FETCH_Q() do { qk = q[q_index]; } while (0)
+    #define BUILD_ROTATION() do {} while(0)
+    #define APPLY_ROTATION() do {} while(0)
+    #define CALL_KERNEL() CALL_IQ(qk,local_values.table)
+
+  // unoriented 1D Beta
+  #elif BETA == 1
+    double F1, F2;
+    #define FETCH_Q() do { qk = q[q_index]; } while (0)
+    #define BUILD_ROTATION() do {} while(0)
+    #define APPLY_ROTATION() do {} while(0)
+    #define CALL_KERNEL() CALL_IQ(qk,F1,F2,local_values.table)
+  #endif
 
 #elif defined(CALL_IQ_A)
@@ -647 +669 @@
       pd_norm += weight * CALL_VOLUME(local_values.table);
       BUILD_ROTATION();
-
+#if BETA
+    if (weight > cutoff) {
+      weight_norm += weight;}
+#endif
 #ifndef USE_OPENCL
       // DLL needs to explicitly loop over the q values.
@@ -693 +718 @@
           }
         #else  // !MAGNETIC
-          const double scattering = CALL_KERNEL();
+          #if defined(CALL_IQ) && BETA == 1
+            CALL_KERNEL();
+            const double scatteringF1 = F1;
+            const double scatteringF2 = F2;
+          #else
+            const double scattering = CALL_KERNEL();
+          #endif
         #endif // !MAGNETIC
 //printf("q_index:%d %g %g %g %g\n", q_index, scattering, weight0);
 
         #ifdef USE_OPENCL
-          this_result += weight * scattering;
+          #if defined(CALL_IQ)&& BETA == 1
+             this_result += weight * scatteringF2;
+             this_beta_result += weight * scatteringF1;
+            #else
+              this_result += weight * scattering;
+          #endif
         #else // !USE_OPENCL
-          result[q_index] += weight * scattering;
+          #if defined(CALL_IQ)&& BETA == 1
+            result[q_index] += weight * scatteringF2;
+            beta_result[q_index] += weight * scatteringF1;
+            #endif
+            #else
+            result[q_index] += weight * scattering;
+          #endif
         #endif // !USE_OPENCL
       }
     }
   }
-
 // close nested loops
 ++step;
@@ -728 +769 @@
   result[q_index] = this_result;
   if (q_index == 0) result[nq] = pd_norm;
+  #if BETA
+  beta_result[q_index] = this_beta_result;
+  #endif
+  if (q_index == 0) result[nq] = pd_norm;
+
 //if (q_index == 0) printf("res: %g/%g\n", result[0], pd_norm);
 #else // !USE_OPENCL
   result[nq] = pd_norm;
+  #if BETA
+  result[2*nq+1] = weight_norm;
+  #endif
 //printf("res: %g/%g\n", result[0], pd_norm);
 #endif // !USE_OPENCL