Minor debugging

spyro/io/model_parameters.py

@@ -551,7 +551,7 @@ def _sanitize_comm(self, comm):
             self.shot_ids_per_propagation = dictionary["parallelism"]["shot_ids_per_propagation"]
         elif self.parallelism_type == "automatic":
             available_cores = COMM_WORLD.size
-            self.shot_ids_per_propagation = [[i] for i in range(0, available_cores)]
+            self.shot_ids_per_propagation = [[i] for i in range(0, self.number_of_sources)]
         elif self.parallelism_type == "spatial":
             self.shot_ids_per_propagation = [[i] for i in range(0, self.number_of_sources)]
 
@@ -633,9 +633,12 @@ def _sanitize_optimization_and_velocity(self):
     def _sanitize_optimization_and_velocity_for_fwi(self):
         self._sanitize_optimization_and_velocity_without_fwi()
         dictionary = self.input_dictionary
-        self.initial_velocity_model_file = dictionary["inversion"][
-            "initial_guess_model_file"
-        ]
+        try:
+            self.initial_velocity_model_file = dictionary["inversion"][
+                "initial_guess_model_file"
+            ]
+        except:
+            self.initial_velocity_model_file = None
         self.fwi_output_folder = "fwi/"
         self.control_output_file = self.fwi_output_folder + "control"
         self.gradient_output_file = self.fwi_output_folder + "gradient"
@@ -666,7 +669,6 @@ def _sanitize_optimization_and_velocity_for_fwi(self):
         if "shot_record_file" in dictionary["inversion"]:
             if dictionary["inversion"]["shot_record_file"] is not None:
                 self.real_shot_record = np.load(dictionary["inversion"]["shot_record_file"])
-
     def _sanitize_optimization_and_velocity_without_fwi(self):
         dictionary = self.input_dictionary
         if "synthetic_data" in dictionary:

test_camembert_fwi.py

@@ -1,22 +1,42 @@
-# from mpi4py.MPI import COMM_WORLD
-# import debugpy
-# debugpy.listen(3000 + COMM_WORLD.rank)
-# debugpy.wait_for_client()
 import spyro
 import firedrake as fire
 import numpy as np
 import matplotlib.pyplot as plt
 import warnings
+import sys
 warnings.filterwarnings("ignore")
 
 
+frequency = float(sys.argv[1])
+degree = int(sys.argv[2])
+
+def cells_per_wavelength(degree):
+    cell_per_wavelength_dictionary = {
+        'kmv2tri': 7.20,
+        'kmv3tri': 3.97,
+        'kmv4tri': 2.67,
+        'kmv5tri': 2.03,
+        'kmv6tri': 1.5,
+        'kmv2tet': 6.12,
+        'kmv3tet': 3.72,
+    }
+
+    cell_type = 'tri'
+
+    key = 'kmv'+str(degree)+cell_type
+
+    return cell_per_wavelength_dictionary.get(key)
+
+cpw = cells_per_wavelength(degree)
 final_time = 0.9
+# dx = 2.5/(frequency*cpw)
+dx = 0.02
 
 dictionary = {}
 dictionary["options"] = {
     "cell_type": "T",  # simplexes such as triangles or tetrahedra (T) or quadrilaterals (Q)
     "variant": "lumped",  # lumped, equispaced or DG, default is lumped
-    "degree": 4,  # p order
+    "degree": degree,  # p order
     "dimension": 2,  # dimension
 }
 dictionary["parallelism"] = {
@@ -33,15 +53,15 @@
     "source_type": "ricker",
     "source_locations": spyro.create_transect((-0.55, 0.7), (-0.55, 1.3), 6),
     # "source_locations": [(-1.1, 1.5)],
-    "frequency": 5.0,
+    "frequency": frequency,
     "delay": 0.2,
     "delay_type": "time",
     "receiver_locations": spyro.create_transect((-1.45, 0.7), (-1.45, 1.3), 200),
 }
 dictionary["time_axis"] = {
     "initial_time": 0.0,  # Initial time for event
     "final_time": final_time,  # Final time for event
-    "dt": 0.001,  # timestep size
+    "dt": 0.00025,  # timestep size
     "amplitude": 1,  # the Ricker has an amplitude of 1.
     "output_frequency": 100,  # how frequently to output solution to pvds - Perguntar Daiane ''post_processing_frequnecy'
     "gradient_sampling_frequency": 1,  # how frequently to save solution to RAM    - Perguntar Daiane 'gradient_sampling_frequency'
@@ -67,7 +87,7 @@ def test_real_shot_record_generation_parallel():
 
     fwi = spyro.FullWaveformInversion(dictionary=dictionary)
 
-    fwi.set_real_mesh(mesh_parameters={"dx": 0.1})
+    fwi.set_real_mesh(mesh_parameters={"dx": dx})
     center_z = -1.0
     center_x = 1.0
     mesh_z = fwi.mesh_z
@@ -80,10 +100,11 @@ def test_real_shot_record_generation_parallel():
 
 def test_realistic_fwi():
     dictionary["inversion"] = {
+        "perform_fwi": True,
         "real_shot_record_files": "shots/shot_record_",
     }
     fwi = spyro.FullWaveformInversion(dictionary=dictionary)
-    fwi.set_guess_mesh(mesh_parameters={"dx": 0.1})
+    fwi.set_guess_mesh(mesh_parameters={"dx": dx})
     fwi.set_guess_velocity_model(constant=2.5)
     mask_boundaries = {
         "z_min": -1.3,
@@ -92,9 +113,9 @@ def test_realistic_fwi():
         "x_max": 1.3,
     }
     fwi.set_gradient_mask(boundaries=mask_boundaries)
-    fwi.run_fwi(vmin=2.5, vmax=3.0, maxiter=5)
+    fwi.run_fwi_rol(vmin=2.5, vmax=3.0, maxiter=30)
 
 
 if __name__ == "__main__":
-    test_real_shot_record_generation_parallel()
-    # test_realistic_fwi()
+    # test_real_shot_record_generation_parallel()
+    test_realistic_fwi()