Feature update interface

pywit/materials.py

@@ -1,4 +1,5 @@
-from pywit.interface import Layer
+from IW2D import IW2DLayer
+from pywit.interface import standard_layer
 from pywit.utils import round_sigfigs
 
 from pathlib import Path
@@ -7,7 +8,7 @@
 from typing import Callable, Tuple
 
 
-def layer_from_dict(thickness: float, material_dict: dict) -> Layer:
+def layer_from_dict(thickness: float, material_dict: dict) -> IW2DLayer:
     """
     Define a layer from a dictionary containing the materials properties.
 
@@ -31,7 +32,7 @@ def layer_from_dict(thickness: float, material_dict: dict) -> Layer:
     assert not any(missing_properties_list), '{} missing from the input dictionary'.format(
             ", ".join(required_material_properties[np.asarray(missing_properties_list)]))
 
-    return Layer(thickness=thickness,
+    return standard_layer(thickness=thickness,
                  dc_resistivity=material_dict['dc_resistivity'],
                  resistivity_relaxation_time=material_dict['resistivity_relaxation_time'],
                  re_dielectric_constant=material_dict['re_dielectric_constant'],
@@ -40,7 +41,7 @@ def layer_from_dict(thickness: float, material_dict: dict) -> Layer:
 
 
 def layer_from_json_material_library(thickness: float, material_key: str,
-                                     library_path: Path = Path(__file__).parent.joinpath('materials.json')) -> Layer:
+                                     library_path: Path = Path(__file__).parent.joinpath('materials.json')) -> IW2DLayer:
     """
     Define a layer using the materials.json library of materials properties.
 
@@ -132,7 +133,7 @@ def magnetoresistance_Kohler(B_times_Sratio: float, P: Tuple[float, ...]) -> flo
         return 10.**np.polyval(P, np.log10(B_times_Sratio))
 
 
-def copper_at_temperature(thickness: float, T: float = 300, RRR: float = 70, B: float = 0) -> Layer:
+def copper_at_temperature(thickness: float, T: float = 300, RRR: float = 70, B: float = 0) -> IW2DLayer:
     """
     Define a layer of pure copper material at any temperature, any B field and any RRR.
     We use a magnetoresistance law fitted from the UPPER curve of the plot in NIST, "Properties of copper and copper
@@ -178,15 +179,15 @@ def copper_at_temperature(thickness: float, T: float = 300, RRR: float = 70, B:
     n = 6.022e23*Z*rho_m*1e6 / A
     tauAC = round_sigfigs(me / (n*dc_resistivity*e**2),3)  # relaxation time (s) (3 significant digits)
 
-    return Layer(thickness=thickness,
+    return standard_layer(thickness=thickness,
                  dc_resistivity=dc_resistivity,
                  resistivity_relaxation_time=tauAC,
                  re_dielectric_constant=1,
                  magnetic_susceptibility=0,
                  permeability_relaxation_frequency=np.inf)
 
 
-def tungsten_at_temperature(thickness: float, T: float = 300, RRR: float = 70, B: float = 0) -> Layer:
+def tungsten_at_temperature(thickness: float, T: float = 300, RRR: float = 70, B: float = 0) -> IW2DLayer:
     """
     Define a layer of tungsten at any temperature, any B field and any RRR.
     The resistivity vs. temperature and RRR is found from Hust & Lankford (see above).
@@ -232,7 +233,7 @@ def tungsten_at_temperature(thickness: float, T: float = 300, RRR: float = 70, B
     n = 6.022e23 * Z * rhom * 1e6 / A
     tauAC = round_sigfigs(me/(n*dc_resistivity*e**2),3) # relaxation time (s) (3 significant digits)
 
-    return Layer(thickness=thickness,
+    return standard_layer(thickness=thickness,
                  dc_resistivity=dc_resistivity,
                  resistivity_relaxation_time=tauAC,
                  re_dielectric_constant=1,

pywit/utilities.py

@@ -1,6 +1,6 @@
 from pywit.component import Component
 from pywit.element import Element
-from pywit.interface import Layer, FlatIW2DInput, RoundIW2DInput
+from IW2D import IW2DLayer, Eps1FromResistivity, Mu1FromSusceptibility, FlatIW2DInput, RoundIW2DInput
 from pywit.interface import component_names
 
 from yaml import load, SafeLoader
@@ -207,7 +207,7 @@ def create_many_resonators_element(length: float, beta_x: float, beta_y: float,
 
 
 def create_classic_thick_wall_component(plane: str, exponents: Tuple[int, int, int, int],
-                                        layer: Layer, radius: float) -> Component:
+                                        layer: IW2DLayer, radius: float) -> Component:
     """
     Creates a single component object modeling a resistive wall impedance/wake,
     based on the "classic thick wall formula" (see e.g. A. W. Chao, chap. 2 in
@@ -220,10 +220,16 @@ def create_classic_thick_wall_component(plane: str, exponents: Tuple[int, int, i
     :param radius: the chamber radius in m
     :return: A component object
     """
+
+    if not isinstance(layer.eps1, Eps1FromResistivity):
+        raise ValueError(f"Input files can only be generated using IW2D.Eps1FromResistivity as eps1 attribute of IW2D.IW2DLayer")
+
+    if not isinstance(layer.mu1, Mu1FromSusceptibility):
+        raise ValueError(f"Input files can only be generated using IW2D.Mu1FromSusceptibility as mu1 attribute of IW2D.IW2DLayer")
 
     # Material properties required for the skin depth computation are derived from the input Layer attributes
-    material_resistivity = layer.dc_resistivity
-    material_relative_permeability = 1. + layer.magnetic_susceptibility
+    material_resistivity = layer.eps1.dc_resistivity
+    material_relative_permeability = 1. + layer.mu1.magnetic_susceptibility
     material_permeability = material_relative_permeability * mu0
 
     # Create the skin depth as a function offrequency and layer properties
@@ -250,7 +256,7 @@ def create_classic_thick_wall_component(plane: str, exponents: Tuple[int, int, i
 
 
 def _zlong_round_single_layer_approx(frequencies: ArrayLike, gamma: float,
-                                     layer: Layer, radius: float, length: float) -> ArrayLike:
+                                     layer: IW2DLayer, radius: float, length: float) -> ArrayLike:
     """
     Function to compute the longitudinal resistive-wall impedance from
     the single-layer, approximated formula for a cylindrical structure,
@@ -271,10 +277,8 @@ def _zlong_round_single_layer_approx(frequencies: ArrayLike, gamma: float,
     omega = 2*pi*frequencies
     k = omega/(beta*c_light)
 
-    rho = layer.dc_resistivity
-    tau = layer.resistivity_relaxation_time
-    mu1 = 1.+layer.magnetic_susceptibility
-    eps1 = 1. - 1j/(eps0*rho*omega*(1.+1j*omega*tau))
+    mu1 = layer.mu1(frequencies)
+    eps1 = layer.eps1(frequencies)
     nu = k*sqrt(1.-beta**2*eps1*mu1)
 
     coef_long = 1j*omega*mu0*length/(2.*pi*beta**2*gamma**2)
@@ -289,7 +293,7 @@ def _zlong_round_single_layer_approx(frequencies: ArrayLike, gamma: float,
 
 
 def _zdip_round_single_layer_approx(frequencies: ArrayLike, gamma: float,
-                                    layer: Layer, radius: float, length: float) -> ArrayLike:
+                                    layer: IW2DLayer, radius: float, length: float) -> ArrayLike:
     """
     Function to compute the transverse dipolar resistive-wall impedance from
     the single-layer, approximated formula for a cylindrical structure,
@@ -310,10 +314,8 @@ def _zdip_round_single_layer_approx(frequencies: ArrayLike, gamma: float,
     omega = 2*pi*frequencies
     k = omega/(beta*c_light)
 
-    rho = layer.dc_resistivity
-    tau = layer.resistivity_relaxation_time
-    mu1 = 1.+layer.magnetic_susceptibility
-    eps1 = 1. - 1j/(eps0*rho*omega*(1.+1j*omega*tau))
+    mu1 = layer.mu1(frequencies)
+    eps1 = layer.eps1(frequencies)
     nu = k*sqrt(1.-beta**2*eps1*mu1)
 
     coef_dip = 1j*k**2*Z0*length/(4.*pi*beta*gamma**4)
@@ -377,11 +379,11 @@ def create_resistive_wall_single_layer_approx_component(plane: str, exponents: T
         if len(input_data.layers) > 1:
             raise NotImplementedError("Input data can have only one layer")
         layer = input_data.layers[0]
-        yok_long = input_data.yokoya_factors[0]
-        yok_dipx = input_data.yokoya_factors[1]
-        yok_dipy = input_data.yokoya_factors[2]
-        yok_quax = input_data.yokoya_factors[3]
-        yok_quay = input_data.yokoya_factors[4]
+        yok_long = input_data.yokoya_zlong
+        yok_dipx = input_data.yokoya_zxdip
+        yok_dipy = input_data.yokoya_zydip
+        yok_quax = input_data.yokoya_zxquad
+        yok_quay = input_data.yokoya_zyquad
     else:
         raise NotImplementedError("Input of type neither FlatIW2DInput nor RoundIW2DInput cannot be handled")
 

setup.py

@@ -41,6 +41,7 @@ def run(self):
         "sortednp",
         "pyyaml",
         "joblib",
+        "IW2D",
     ],
     "test":[
         "pytest",

test/test_interface.py

@@ -2,11 +2,12 @@
 
 import pytest
 
-from pywit.interface import import_data_iw2d, create_component_from_data, Sampling
+from pywit.interface import component_recipes_from_legacy_iw2d_files, create_component_from_data
 from pywit.interface import check_valid_working_directory, get_component_name
 from pywit.interface import check_already_computed, get_iw2d_config_value, RoundIW2DInput, FlatIW2DInput, add_iw2d_input_to_database
 from pywit.parameters import *
 from pywit.materials import layer_from_json_material_library
+from IW2D import RoundIW2DInput, FlatIW2DInput, InputFileFreqParams, InputFileWakeParams
 
 from pathlib import Path
 import glob
@@ -41,7 +42,7 @@ def test_get_component_name_raise(is_impedance, plane, exponents):
 
 def test_duplicate_component_iw2d_import():
     with raises(AssertionError) as error_message:
-        import_data_iw2d(directory=Path("test/test_data/iw2d/duplicate_components").resolve(),
+        component_recipes_from_legacy_iw2d_files(directory=Path("test/test_data/iw2d/duplicate_components").resolve(),
                          common_string="WLHC_2layersup_0layersdown6.50mm")
 
     assert error_message.value.args[0] in ["The wake files 'WlongWLHC_2layersup_0layersdown6.50mm.dat' and "
@@ -55,7 +56,7 @@ def test_duplicate_component_iw2d_import():
 
 def test_no_matching_filename_iw2d_import():
     with raises(AssertionError) as error_message:
-        import_data_iw2d(directory=Path("test/test_data/iw2d/valid_directory").resolve(),
+        component_recipes_from_legacy_iw2d_files(directory=Path("test/test_data/iw2d/valid_directory").resolve(),
                          common_string="this_string_matches_no_file")
 
     expected_error_message = f"No files in " \
@@ -69,13 +70,11 @@ def test_valid_iw2d_component_import():
     # Normally, the relativistic gamma would be an attribute of a required IW2DInput object, but here it has been
     # hard-coded instead
     relativstic_gamma = 479.605064966
-    recipes = import_data_iw2d(directory=Path("test/test_data/iw2d/valid_directory").resolve(),
+    recipes = component_recipes_from_legacy_iw2d_files(directory=Path("test/test_data/iw2d/valid_directory").resolve(),
                                common_string="precise")
     for recipe in recipes:
         component = create_component_from_data(*recipe, relativistic_gamma=relativstic_gamma)
-        data = recipe[-1]
-        x = data[:, 0]
-        y_actual = data[:, 1] + (1j * data[:, 2] if data.shape[1] == 3 else 0)
+        x, y_actual = recipe[-2:]
         if recipe[0]:
             np.testing.assert_allclose(y_actual, component.impedance(x), rtol=REL_TOL, atol=ABS_TOL)
         else:
@@ -84,20 +83,29 @@ def test_valid_iw2d_component_import():
 
 @pytest.fixture
 def iw2d_input(request):
-    f_params = Sampling(start=1, stop=1e9, scan_type=0, added=(1e2,))
-    z_params = Sampling(start=1e-9, stop=1, scan_type=0, added=(1e-6,))
+    f_params = InputFileFreqParams(use_log_sampling=True, use_lin_sampling=False, log_fmin=1, log_fmax=1e9, log_f_per_decade=1, added_f=(1e2,))
+    wake_params = InputFileWakeParams(long_error_weight=1, wake_abs_tolerance=1, freq_lin_bisect=1e9,
+                                      use_log_sampling=True, use_lin_sampling=False, log_zmin=1e-9, log_zmax=1, log_z_per_decade=1, added_z=(1e-6,))
     layers_tung = (layer_from_json_material_library(thickness=np.inf, material_key='W'),)
 
     if request.param['chamber_type'] == 'round':
-        return RoundIW2DInput(machine='test', length=1, relativistic_gamma=7000,
-                              calculate_wake=request.param['wake_computation'], f_params=f_params, comment='test',
-                              layers=layers_tung, inner_layer_radius=5e-2, yokoya_factors=(1, 1, 1, 1, 1), z_params=z_params)
+        input_object = RoundIW2DInput(length=1, relativistic_gamma=7000,
+                                      layers=layers_tung, inner_layer_radius=5e-2,
+                                      yokoya_zlong = 1, yokoya_zxdip = 1,
+                                      yokoya_zydip = 1, yokoya_zxquad = 0,
+                                      yokoya_zyquad = 0)
 
     if request.param['chamber_type'] == 'flat':
-        return FlatIW2DInput(machine='test', length=1, relativistic_gamma=7000,
-                             calculate_wake=request.param['wake_computation'], f_params=f_params, comment='test',
-                             top_bottom_symmetry=True, top_layers=layers_tung, top_half_gap=5e-2, z_params=z_params)
-
+        input_object = FlatIW2DInput(length=1, relativistic_gamma=7000,
+                             top_bottom_symmetry=True, top_layers=layers_tung, top_half_gap=5e-2)
+
+    if request.param['wake_computation'] == True:
+        additional_params = wake_params
+    else:
+        additional_params = f_params
+
+    return (input_object, additional_params)
+
 
 def _remove_non_empty_directory(directory_path: Path):
     if not os.path.exists(directory_path):
@@ -116,11 +124,12 @@ def _remove_non_empty_directory(directory_path: Path):
                           ({'chamber_type': 'flat', 'wake_computation': True}, ['Zlong', 'Zxdip', 'Zydip', 'Zxquad', 'Zyquad', 'Zycst', 'Wlong', 'Wxdip', 'Wydip', 'Wxquad', 'Wyquad', 'Wycst'])]
 @pytest.mark.parametrize("iw2d_input, components_to_test", list_of_inputs_to_test, indirect=["iw2d_input"])
 def test_check_already_computed(iw2d_input, components_to_test):
+    iw2d_input, additional_input_params = iw2d_input
     name = 'test_hash'
 
     # create the expected directories for the dummy input
     projects_path = Path(get_iw2d_config_value('project_directory'))
-    input_hash = sha256(iw2d_input.__str__().encode()).hexdigest()
+    input_hash = iw2d_input.input_file_hash(additional_input_params)
     directory_level_1 = projects_path.joinpath(input_hash[0:2])
     directory_level_2 = directory_level_1.joinpath(input_hash[2:4])
     working_directory = directory_level_2.joinpath(input_hash[4:])
@@ -143,18 +152,18 @@ def test_check_already_computed(iw2d_input, components_to_test):
             f.write(dummy_string)
 
     # check that the input is detected in the hashmap
-    already_computed, input_hash, working_directory = check_already_computed(iw2d_input, name)
+    already_computed, input_hash, working_directory = check_already_computed(iw2d_input, additional_input_params, name)
     assert already_computed
 
     # now we remove one component and verify that check_already_computed gives false
     os.remove(f'{working_directory}/Zlong_test.txt')
-    already_computed, input_hash, working_directory = check_already_computed(iw2d_input, name)
+    already_computed, input_hash, working_directory = check_already_computed(iw2d_input, additional_input_params, name)
     assert not already_computed
 
     # now we remove the folder and check that check_already_computed gives false
     _remove_non_empty_directory(working_directory)
 
-    already_computed, input_hash, working_directory = check_already_computed(iw2d_input, name)
+    already_computed, input_hash, working_directory = check_already_computed(iw2d_input, additional_input_params, name)
 
     assert not already_computed
 
@@ -164,13 +173,14 @@ def test_check_already_computed(iw2d_input, components_to_test):
 
 @pytest.mark.parametrize("iw2d_input", [{'chamber_type': 'round', 'wake_computation': False}], indirect=["iw2d_input"])
 def test_add_iw2d_input_to_database(iw2d_input):
+    iw2d_input, additional_input_params = iw2d_input
     projects_path = Path(get_iw2d_config_value('project_directory'))
-    input_hash = sha256(iw2d_input.__str__().encode()).hexdigest()
+    input_hash = iw2d_input.input_file_hash(additional_input_params)
     directory_level_1 = projects_path.joinpath(input_hash[0:2])
     directory_level_2 = directory_level_1.joinpath(input_hash[2:4])
     working_directory = directory_level_2.joinpath(input_hash[4:])
 
-    add_iw2d_input_to_database(iw2d_input, input_hash, working_directory)
+    add_iw2d_input_to_database(iw2d_input, additional_input_params, input_hash, working_directory)
 
     assert os.path.exists(f"{working_directory}/input.txt")
 

test/test_materials.py

@@ -40,19 +40,19 @@ def test_layer_from_dict_two_missing_keys():
 def test_copper(T,B,rho):
     RRR = 70.
     copper_layer = materials.copper_at_temperature(1., T, RRR, B)
-    testing.assert_allclose(copper_layer.dc_resistivity, rho, atol=1e-11)
+    testing.assert_allclose(copper_layer.eps1.dc_resistivity, rho, atol=1e-11)
 
 
 @mark.parametrize("material_key, material_function, T, RRR",
                   [ ["Cu", 'copper_at_temperature', 293, 70],
                     ["W",  'tungsten_at_temperature', 300, 70],
                   ])
 @mark.parametrize("material_property, tolerance",
-                  [ ['dc_resistivity', 2e-9],
-                    ['resistivity_relaxation_time', 3e-15], 
-                    ['re_dielectric_constant', 0],
-                    ['magnetic_susceptibility', 0],
-                    ['permeability_relaxation_frequency', 0],
+                  [ ['eps1.dc_resistivity', 2e-9],
+                    ['eps1.resistivity_relaxation_time', 3e-15], 
+                    ['eps1.re_dielectric_constant', 0],
+                    ['mu1.magnetic_susceptibility', 0],
+                    ['mu1.permeability_relaxation_frequency', 0],
                   ])
 def test_materials_at_temperature(material_key, material_function, T,
                                   RRR, material_property, tolerance):
@@ -61,7 +61,15 @@ def test_materials_at_temperature(material_key, material_function, T,
     layer_from_library = materials.layer_from_json_material_library(1.,material_key)
 
     assert layer_at_temperature.thickness == layer_from_library.thickness == 1.
-    testing.assert_allclose(getattr(layer_at_temperature,material_property),
-                            getattr(layer_from_library,material_property),
-                            atol=tolerance)
+
+    eps_or_mu, property_name = material_property.split(".")
+    if eps_or_mu == "eps1":
+        testing.assert_allclose(getattr(layer_at_temperature.eps1,property_name),
+                                getattr(layer_from_library.eps1,property_name),
+                                atol=tolerance)
+    else:
+      testing.assert_allclose(getattr(layer_at_temperature.mu1,property_name),
+                                getattr(layer_from_library.mu1,property_name),
+                                atol=tolerance)
+