Set up the Gaussians for the Li2 example

LLNL · Aug 8, 2024 · 4742344 · 4742344
1 parent b42ab46
commit 4742344
Showing 1 changed file with 11 additions and 16 deletions.
diff --git a/examples/prom/elliptic_eigenproblem_global_rom.cpp b/examples/prom/elliptic_eigenproblem_global_rom.cpp
@@ -852,6 +852,7 @@ double Potential(const Vector &x)
 {
     Vector center(x.Size());
     double rho = 0.0;
+    double scaling_factor = pow((bb_max[0] - bb_min[0]) / 18.0, 2.0);
     switch (problem)
     {
     case 1:
@@ -868,26 +869,20 @@ double Potential(const Vector &x)
         if (potential_well_switch == 0 || potential_well_switch == 1)
         {
             // add well with first center
-            center(0) = (0.25 + 0.02 * cos(2 * M_PI * pseudo_time)) *
-                        (bb_min[0] + bb_max[0]);
-            center(1) = (0.25 + 0.02 * sin(2 * M_PI * pseudo_time)) *
-                        (bb_min[1] + bb_max[1]);
-            for (int i = 2; i < x.Size(); i++)
-                center(i) = 0.25 * (bb_min[i] + bb_max[i]);
-            rho += amplitude * std::exp(-x.DistanceSquaredTo(center) / pow(4.0 * h_max,
-                                        2.0));
+            center(0) = (23 * bb_min[0] + 13 * bb_max[0]) / 36 - pseudo_time * h_max;
+            for (int i = 1; i < x.Size(); i++)
+                center(i) = 0.0;
+            rho += amplitude * scaling_factor * std::exp(-50 / scaling_factor *
+                    x.DistanceSquaredTo(center));
         }
         if (potential_well_switch == 0 || potential_well_switch == 2)
         {
             // add well with second center
-            center(0) = (0.75 - 0.02 * cos(2 * M_PI * pseudo_time)) *
-                        (bb_min[0] + bb_max[0]);
-            center(1) = (0.75 - 0.02 * sin(2 * M_PI * pseudo_time)) *
-                        (bb_min[1] + bb_max[1]);
-            for (int i = 2; i < x.Size(); i++)
-                center(i) = 0.75 * (bb_min[i] + bb_max[i]);
-            rho += amplitude * std::exp(-x.DistanceSquaredTo(center) / pow(4.0 * h_max,
-                                        2.0));
+            center(0) = (23 * bb_min[0] + 23 * bb_max[0]) / 36 + pseudo_time * h_max;
+            for (int i = 1; i < x.Size(); i++)
+                center(i) = 0.0;
+            rho += amplitude * scaling_factor * std::exp(-50 / scaling_factor *
+                    x.DistanceSquaredTo(center));
         }
         return rho;
     }