Initial reliability analysis

examples/infeas_analysis_cluster.ipynb

@@ -424,7 +424,8 @@
     }
    ],
    "source": [
-    "campaign, samples = analyse.read_campaign(campaign_name=\"example_cluster\")\n",
+    "campaign = analyse.get_campaign(campaign_name=\"example_cluster\")\n",
+    "samples = analyse.get_samples(campaign)\n",
     "samples"
    ]
   },

examples/infeas_analysis_local.ipynb

@@ -411,7 +411,8 @@
     }
    ],
    "source": [
-    "campaign, samples = analyse.read_campaign(campaign_name=\"example_local\")\n",
+    "campaign = analyse.get_campaign(campaign_name=\"example_local\")\n",
+    "samples = analyse.get_samples(campaign)\n",
     "samples"
    ]
   },

infeas/analyse.py

@@ -9,8 +9,8 @@
 from infeas.eval import QOIS, WORK_DIR
 
 
-def read_campaign(campaign_name: str) -> pd.DataFrame:
-    """Read in evaluated campaign and return dataframe.
+def get_campaign(campaign_name: str) -> uq.campaign.Campaign:
+    """Read in evaluated campaign and return campaign object.
 
     Fetches the latest campaign with matching name.
 
@@ -23,8 +23,6 @@ def read_campaign(campaign_name: str) -> pd.DataFrame:
     -------
     uq.campaign.Campaign
         easyVVUQ campaign
-    pd.DataFrame
-        Evaluation results
 
     Raises
     ------
@@ -52,13 +50,27 @@ def read_campaign(campaign_name: str) -> pd.DataFrame:
         db_location=db_location_prefixed, name=campaign_name, work_dir=WORK_DIR
     )
 
+    return campaign
+
+
+def get_samples(campaign: uq.campaign.Campaign) -> pd.DataFrame:
+    """Get samples from campaign object.
+
+    Parameters
+    ----------
+    campaign : uq.campaign.Campaign
+        campaign object
+
+    Returns
+    -------
+    pd.DataFrame
+        evaluated samples
+    """
     samples = campaign.get_collation_result()
-    sample_count = samples.shape[0]
-    print(f"Campaign read in. Number of samples = {sample_count}")
 
-    # Drop strange multi-index of 0
+    # Drop unnecessary multi-index of 0
     samples.columns = samples.columns.droplevel(1)
-    return campaign, samples
+    return samples
 
 
 def describe_qois(results: pd.DataFrame) -> pd.DataFrame:

infeas/decoder.py

@@ -73,27 +73,21 @@ def _process_raw_data(self, raw_data: Dict[str, Any]) -> Dict[str, float]:
             for ineq_constrs_key in ineq_constrs_keys
         }
 
-        # Only want violated constraint values
-        # Coerce feasible inequality constraints (> 0) = 0.0
-        # TODO Not sure if we want to mask non-violated constraint
-        # values at this stage: infeasibile responses only
-        vio_ineq_constrs_dict = {}
-        for key, value in ineq_constrs_dict.items():
-            if value > 0:
-                vio_ineq_constrs_dict[key] = 0.0
-            else:
-                vio_ineq_constrs_dict[key] = value
-
-        # Merge individual eq and ineq violated constraint values
-        responses = responses | eq_constrs_dict | vio_ineq_constrs_dict
+        # Merge individual eq and ineq constraint values
+        responses = responses | eq_constrs_dict | ineq_constrs_dict
 
         # Calculate RMS constraint residuals for violated constraints only
         # Create arrays from constraints dicts
-        eq_constrs = np.array(list(eq_constrs_dict.values()))
-        vio_ineq_constrs = np.array(list(vio_ineq_constrs_dict.values()))
-        vio_constrs = np.concatenate((eq_constrs, vio_ineq_constrs))
-        rms_vio_constr_res = np.sqrt(np.mean(vio_constrs**2))
-        responses["rms_vio_constr_res"] = rms_vio_constr_res
+        # Process satisfied con: c > 0
+        # Everyone else satisfied con: c < 0
+        # Flip sign to agree with standard
+        ineq_constrs = -np.array(list(ineq_constrs_dict.values()))
+
+        # TODO Would be better to override method in notebook for testing
+
+        # Use w metric: most violated inequality constraint
+        w = np.max(ineq_constrs)
+        responses["w"] = w
 
         return responses
 

infeas/mfile_to_in.py

@@ -9,6 +9,7 @@
 from pathlib import Path
 from typing import Optional, Sequence
 import re
+from shutil import copy
 
 
 """Unfortunately, f-values have to be detected in two ways. This is because
@@ -120,34 +121,111 @@
 
 
 def convert(
-    mfile_name: str,
     original_in_name: str,
-    sol_in_name: str,
+    new_in_name: str,
     no_optimisation: Optional[bool] = True,
     n_equalities: Optional[int] = None,
+    mfile_name: Optional[str] = None,
+    remove_f_value_inits: Optional[bool] = False,
 ) -> Path:
-    """Convert mfile to input file, preserving optimisation parameter vector.
+    """Convert input file to a modified new one.
+
+    If input file and n_equalities provided, just convert input to inequalities.
+    If mfile is provided, copy solution vector to new input file.
 
     Parameters
     ----------
-    mfile_name : str
-        name of mfile to convert
     original_in_name : str
         the IN.DAT used to create the MFILE.DAT solution
-    sol_in_name : str
-        name of input file to be created at solution vector
+    new_in_name : str
+        modified input file to be created
     no_optimisation : Optional[bool], optional
         convert IN.DAT to non-optimising run, by default True
     n_equalities : Optional[int], optional
         how many of the constraints to make equalities. The first n_equalities
         constraints will be left as equalities, the rest converted to
         inequalities, by default None
+    mfile_name : Optional[str], optional
+        if mfile supplied, extract solution vector into new input file, by
+        default None
+    remove_f_value_inits : Optional[bool], optional
+        remove f-value initialisations, e.g. `fthresh = 0.8`, by default False
 
     Returns
     -------
     Path
         path to solution input file
     """
+    # Create new input file, optionally overwriting with solution vector
+    if mfile_name:
+        # Extract solution vector from mfile, copy into new input file
+        sol_in_dat_path = overwrite_sol_vector(
+            mfile_name=mfile_name,
+            original_in_name=original_in_name,
+            new_in_name=new_in_name,
+        )
+    else:
+        # Just copy original input to new path
+        sol_in_dat_path = copy(original_in_name, new_in_name)
+
+    # Additional modifications to top of input file
+    top_content = []
+    if no_optimisation:
+        top_content.extend(["* Once through only: no optimisation", "ioptimz  = -2"])
+
+    if n_equalities is not None:
+        # Write neqns count
+        top_content.extend(
+            [
+                "* Define number of equality constraints, and",
+                "* use inequality constraints: corresponding f-values removed",
+                f"neqns = {n_equalities}",
+            ]
+        )
+
+    new_content = [line + "\n" for line in top_content]
+
+    # Write top and main content, removing f-values if necessary
+    with open(sol_in_dat_path, "r") as f:
+        content = f.readlines()
+
+    if n_equalities:
+        # Remove all f-values as opt params
+        main_content = remove_f_values(content, remove_f_value_inits)
+    else:
+        main_content = content
+
+    new_content.extend(main_content)
+
+    with open(sol_in_dat_path, "w") as f:
+        f.writelines(new_content)
+
+    return sol_in_dat_path
+
+
+def overwrite_sol_vector(
+    mfile_name: str, original_in_name: str, new_in_name: str
+) -> Path:
+    """Extract solution vector from mfile, copy into new input file.
+
+    Parameters
+    ----------
+    mfile_name : str
+        output containing solution vector
+    original_in_name : str
+        original input file
+    new_in_name : str
+        input file with initial optimisation parameter vector overwritten with
+        the solution vector
+
+    Returns
+    -------
+    Path
+        path to new solution vector input file
+    """
+    # Pattern for single array element, e.g. fimp(13)
+    array_el_re = re.compile(r"(\w+)\((\d+)\)")
+
     # First, extract the solution (optimisation parameters) from the mfile
     # Create Mfile object from mfile
     mfile_path = Path(mfile_name)
@@ -189,63 +267,53 @@ def convert(
                 continue
 
         # Otherwise update new IN.DAT value with solution value
-        in_dat.add_parameter(var_name, value)
+        matches = array_el_re.match(var_name)
+        if matches is not None:
+            # Found a single array element e.g. fimp(13)
+            array_name = matches.group(1)
+            array_index = matches.group(2)
+            py_array_index = int(array_index) - 1
+
+            in_dat.change_array(
+                array_name=array_name, array_index=py_array_index, array_value=value
+            )
+        else:
+            # Set individual parameter
+            in_dat.add_parameter(var_name, value)
 
-    print(f"Solution f-values ignored = {dropped_sol_f_value_count}")
+    print(f"Optimisation parameter f-values ignored = {dropped_sol_f_value_count}")
 
     # Write out new IN.DAT, with optimisation parameters set to original
     # solution vector
-    sol_in_dat_path = Path(sol_in_name)
+    sol_in_dat_path = Path(new_in_name)
     in_dat.write_in_dat(sol_in_dat_path)
-
-    # Additional modifications to top of input file
-    top_content = []
-    if no_optimisation:
-        top_content.extend(["* Once through only: no optimisation", "ioptimz  = -2"])
-
-    if n_equalities is not None:
-        # Write neqns count
-        top_content.extend(
-            [
-                "* Define number of equality constraints, and",
-                "* use inequality constraints: corresponding f-values removed",
-                f"neqns = {n_equalities}",
-            ]
-        )
-
-    new_content = [line + "\n" for line in top_content]
-
-    # Write top and main content, removing f-values if necessary
-    with open(sol_in_dat_path, "r+") as f:
-        content = f.readlines()
-        f.seek(0, 0)
-
-        if n_equalities:
-            # Remove all f-values
-            main_content = remove_f_values(content)
-        else:
-            main_content = content
-
-        new_content.extend(main_content)
-        f.writelines(new_content)
-
     return sol_in_dat_path
 
 
-def remove_f_values(lines_with_f_values: Sequence[str]) -> list[str]:
+def remove_f_values(
+    lines_with_f_values: Sequence[str], remove_f_value_inits: bool
+) -> list[str]:
     """Remove f-value optimisation parameters from input file lines.
 
+    Removes `ixc = 23 * fthresh` line if opt param 23 is f-value.
+    Optionally removes f-value value initialisations `fthresh = 0.8`.
+
     Parameters
     ----------
     lines_with_f_values : Sequence[str]
         lines from original input file
+    remove_f_value_inits : bool
+        remove f-value initialisations, e.g. `fthresh = 0.8`
 
     Returns
     -------
     list[str]
         lines with f-value optimisation parameters removed
     """
+    # f-value opt params removed (ixc = 23 * fdene)
     f_value_removal_count = 0
+    # f-value initialisations removed (fthresh = 0.8)
+    f_value_init_removal_count = 0
     # Lines that will be included in new IN.DAT
     lines_without_f_values = []
     # Opt params beginning with f that aren't on the ignore list
@@ -255,10 +323,14 @@ def remove_f_values(lines_with_f_values: Sequence[str]) -> list[str]:
     # Don't also search for comment (e.g. ixc = 23 * fdene) as it may not be there:
     # Just look for opt param number
     opt_param_re = re.compile(r"ixc\s*=\s*(\d+)")
+
     # If not already found in the f-value opt param number list, attempt to
     # find any opt params starting with f (depends on "*" comment)
     opt_param_beginning_with_f_re = re.compile(r"ixc\s*=\s*(\d+)\s*\*\s*(f\w+)")
 
+    # Any line starting with f (could be f-value value definition)
+    line_starting_with_f_re = re.compile(r"^(f\w*)")
+
     for line in lines_with_f_values:
         matches = opt_param_re.match(line)
         if matches is not None:
@@ -272,9 +344,6 @@ def remove_f_values(lines_with_f_values: Sequence[str]) -> list[str]:
                 f_value_removal_count += 1
                 continue
 
-        # Non f-value opt param or other line: keep
-        lines_without_f_values.append(line)
-
         # Check if line contains opt param beginning with an f:
         # Not known f-value (would be in numbers list), but suspicious as could be
         matches = opt_param_beginning_with_f_re.match(line)
@@ -285,6 +354,20 @@ def remove_f_values(lines_with_f_values: Sequence[str]) -> list[str]:
                 # Could be an f-value
                 suspect_opt_params_beginning_with_f.append(line)
 
+        # Optionally, remove any f-value initialisations
+        # e.g. `fthresh = 0.8`
+        if remove_f_value_inits:
+            matches = line_starting_with_f_re.match(line)
+            if matches is not None:
+                if matches.group(1) not in STARTS_WITH_F_BUT_NOT_F_VALUE_NAMES:
+                    # Found f-value initialisation: drop the line
+                    f_value_init_removal_count += 1
+                    continue
+
+        # Non f-value opt param, (optionally not) f-value initialisation or other
+        # line: keep
+        lines_without_f_values.append(line)
+
     # Print out potential f-value vars if any
     if len(suspect_opt_params_beginning_with_f) > 0:
         raise ValueError(
@@ -293,4 +376,7 @@ def remove_f_values(lines_with_f_values: Sequence[str]) -> list[str]:
     else:
         print(f"f-values removed as optimisation parameters = {f_value_removal_count}")
 
+    if remove_f_value_inits:
+        print(f"f-value initialisations removed = {f_value_init_removal_count}")
+
     return lines_without_f_values

max_net_elec/max_net_elec_analysis.ipynb

@@ -1247,7 +1247,8 @@
     }
    ],
    "source": [
-    "campaign, results = analyse.read_campaign(campaign_name)\n",
+    "campaign = analyse.get_campaign(campaign_name)\n",
+    "results = analyse.get_samples(campaign)\n",
     "results"
    ]
   },