Update README.md in examples and fix CLI bugs

.pre-commit-config.yaml

@@ -7,8 +7,8 @@ repos:
   - id: trailing-whitespace
   - id: end-of-file-fixer
   - id: check-yaml
-    exclude: ^conda/
   - id: check-added-large-files
+    exclude: ^example/AlN-LDA/
 
 - repo: https://github.com/astral-sh/ruff-pre-commit
   rev: v0.5.6

example/AlN-LDA/README.md

@@ -14,52 +14,56 @@ from all displaced supercell forces.
 
 Perfect and displaced supercells were created by
 
-```
+```bash
 % phono3py --dim 3 3 2 -c POSCAR-unitcell -d
 ```
 
-In the example directory, `FORCES_FC3` is compressed to `FORCES_FC3.lzma`. After
-unzipping `FORCES_FC3.lzma` (e.g., using `tar xvfz` or `tar xvfa`), to obtain
-`fc3.hdf5` and normal `fc2.hdf5`,
-
-```
-% phono3py --sym-fc
-```
-
-Using 13x13x9 sampling mesh, lattice thermal conductivity is calculated by
+In the example directory, `FORCES_FC3` is compressed to `FORCES_FC3.xz`. After
+unzipping `FORCES_FC3.xz` (e.g., using `xz -d`), to obtain `fc3.hdf5` and
+`fc2.hdf5` using symfc (the results without using symfc, i.e., finite difference
+method, are shown at the bottom of this README)
 
-```
-% phono3py --mesh 13 13 9 --fc3 --fc2 --br
+```bash
+% phono3py-load --symfc -v
 ```
 
-`kappa-m13139.hdf5` is written as the result. The lattice thermal conductivity
-is calculated as k_xx=228.2 and k_zz=224.1 W/m-K at 300 K.
+Lattice thermal conductivity is calculated by
 
-With `--nac` option, non-analytical term correction is applied reading the Born
-effective charges and dielectric constant from `BORN` file:
-
-```
-% phono3py --mesh 13 13 9 --fc3 --fc2 --br --nac
+```bash
+% phono3py-load --mesh 40 --br --ts 300
 ```
 
-This changes thermal conductivity at 300 K to k_xx=235.7 and k_zz=219.1. The
-shape of phonon band structure is important to fullfil energy and momentum
-conservations.
+`kappa-m15158.hdf5` is written as the result. Parameters for non-analytical term
+correction (NAC) is automatically read from those stored in `phono3py_disp.yaml` or
+`BORN` file. The lattice thermal conductivity is calculated as k_xx=242.8 and
+k_zz=226.5 W/m-K at 300 K. Without NAC, k_xx=233.6 and k_zz=222.2.
 
-Use of larger supercell of fc2 may change the shape of phonon band structure. To
-see it, first regenerate `phono3py_disp.yaml` with `--dim-fc2` option,
+Use of larger supercell for fc2 may change the shape of phonon band structure.
+To see it, first regenerate `phono3py_disp.yaml` with `--dim-fc2` option,
 
-```
+```bash
 % phono3py --dim 3 3 2 --dim-fc2 5 5 3 -c POSCAR-unitcell -d
 ```
 
 Then re-create force constants and calculate thermal conductivity,
 
+```bash
+% phono3py-load --symfc -v
+% phono3py-load --br --mesh=40 --ts 300
 ```
-% phono3py --sym-fc
-% phono3py --mesh="13 13 9" --fc3 --fc2 --br --nac
+
+If `phono3py_disp.yaml` is renamed to `phono3py_disp_dimfc2.yaml`, it can be
+specified at the first argument of `phono3py-load` command:
+
+```bash
+% phono3py-load phono3py_disp_dimfc2.yaml --symfc -v
+% phono3py-load phono3py_disp_dimfc2.yaml --br --mesh=40 --ts 300
 ```
 
-k_xx=236.0 and k_zz=222.2 are obtained. In the case of this example, we can see
+k_xx=240.2 and k_zz=230.1 are obtained. In the case of this example, we can see
 that the larger fc2 supercell contributes little, which means that the 3x3x2
 supercell was good enough to obtain a good shape of phonon band structure.
+
+Using the finite difference method implemented in phono3py, lattice thermal
+conductivities are obtained as k_xx=251.2 and k_zz=233,4 without using the large
+fc2 supercell and k_xx=249.4 k_zz=236.9 using the large fc2 supercell.

example/NaCl-alm/README.md

@@ -29,4 +29,4 @@ lattice thermal conductivity at 300 K is calculated by
 % phono3py-load phono3py_params_NaCl222.yaml.xz --mesh 50 --ts 300 --br
 ```
 
-The result is ~7.2 W/m-K at 300 K.
+The result is ~7.3 W/m-K at 300 K.

example/Si-LDA/README.md

@@ -13,18 +13,18 @@ conventional unit cell to the primitive cell.
 To create `fc3.hdf5` and `fc2.hdf5`,
 
 ```
-% phono3py --sym-fc
+% phono3py-load
 ```
 
 Using 11x11x11 sampling mesh, lattice thermal conductivity is calculated by
 
 ```
-% phono3py --mesh 11 11 11 --fc3 --fc2 --br
+% phono3py-load --mesh 11 11 11 --br --ts 300
 ```
 
 `kappa-m111111.hdf5` is written as the result. The lattice thermal conductivity
-is calculated as 112.5 W/m-K at 300 K. This becomes, with 19x19x19 sampling
-mesh, 127.0 W/m-K.
+is calculated as 112.4 W/m-K at 300 K. This becomes, with 19x19x19 sampling
+mesh, 128.2 W/m-K.
 
 Accumulated lattice thermal conductivity is calculated using `phono3py-kaccum`
 script.

example/Si-PBE/README.md

@@ -13,18 +13,18 @@ conventional unit cell to the primitive cell.
 To create `fc3.hdf5` and `fc2.hdf5`,
 
 ```
-% phono3py --sym-fc
+% phono3py-load
 ```
 
 Using 11x11x11 sampling mesh, lattice thermal conductivity is calculated by
 
 ```
-% phono3py --mesh 11 11 11 --fc3 --fc2 --br
+% phono3py-load --mesh 11 11 11 --fc3 --fc2 --br
 ```
 
 `kappa-m111111.hdf5` is written as the result. The lattice thermal conductivity
-is calculated as 119.3 W/m-K at 300 K. This becomes, with 19x19x19 sampling
-mesh, 132.4 W/m-K.
+is calculated as 119.5 W/m-K at 300 K. This becomes, with 19x19x19 sampling
+mesh, 130.1 W/m-K.
 
 Accumulated lattice thermal conductivity is calculated using `phono3py-kaccum`
 script.

example/Si-PBEsol/README.md

@@ -13,18 +13,18 @@ the conventional unit cell to the primitive cell.
 To create `fc3.hdf5` and `fc2.hdf5`,
 
 ```
-% phono3py --sym-fc
+% phono3py-load
 ```
 
 Using 11x11x11 sampling mesh, lattice thermal conductivity is calculated by
 
 ```
-% phono3py --mesh 11 11 11 --fc3 --fc2 --br
+% phono3py-load --mesh 11 11 11 --br
 ```
 
 `kappa-m111111.hdf5` is written as the result. The lattice thermal conductivity
-is calculated as 109.0 W/m-K at 300 K. This becomes, with 19x19x19 sampling
-mesh, 123.2 W/m-K.
+is calculated as 109.1 W/m-K at 300 K. This becomes, with 19x19x19 sampling
+mesh, 124.4 W/m-K.
 
 Accumulated lattice thermal conductivity is calculated using `phono3py-kaccum`
 script.
@@ -80,11 +80,11 @@ In this example directory, the dataset is renamed to
 at `vasprun_xml_fc2/disp-fc2-00001/vasprun.xml`. `FORCES_FC2` is generated by
 
 ```
-% phono3py phono3py_disp_dimfc2.yaml --cf2 vasprun_xml_fc2/disp-fc2-00001/vasprun.xml
+% phono3py -c phono3py_disp_dimfc2.yaml --cf2 vasprun_xml_fc2/disp-fc2-00001/vasprun.xml
 ```
 
 A lattice thermal conductivity calculation is performed by, e.g.,
 ```
 % phono3py-load phono3py_disp_dimfc2.yaml --mesh 11 11 11 --br --ts 300
 ```
-The result is 108.0 W/m-K, and with the 19x19x19 mesh, 125.4 W/m-K.
+The result is 107.9 W/m-K, and with the 19x19x19 mesh, 125.4 W/m-K.

example/Si-QE/README.md

@@ -13,18 +13,18 @@ the conventional unit cell to the primitive cell.
 To create `fc3.hdf5` and `fc2.hdf5`,
 
 ```
-% phono3py --sym-fc
+% phono3py-load
 ```
 
 Using 11x11x11 sampling mesh, lattice thermal conductivity is calculated by
 
 ```
-% phono3py --mesh 11 11 11 --fc3 --fc2 --br
+% phono3py-load --mesh 11 11 11 --br --ts 300
 ```
 
 `kappa-m111111.hdf5` is written as the result. The lattice thermal conductivity
-is calculated as 118.9 W/m-K at 300 K. This becomes, with 19x19x19 sampling
-mesh, 129.9 W/m-K.
+is calculated as 121.3 W/m-K at 300 K. This becomes, with 19x19x19 sampling
+mesh, 135.4 W/m-K.
 
 The .out files for supercells are found in `supercell_out.tar.lzma`. If phono3py
 is properly installed, the following command should work.

phono3py/api_phono3py.py

@@ -1487,6 +1487,9 @@ def produce_fc3(
             fc2 = None
         if fc_calculator_options is not None and "|" in fc_calculator_options:
             fc2 = None
+        # fc2 should not be set if phonon_supercell_matrix is available.
+        if self._phonon_supercell_matrix is not None:
+            fc2 = None
 
         # Normally self._fc2 is overwritten in produce_fc2
         if self._fc2 is None:
@@ -1527,6 +1530,9 @@ def produce_fc2(
         else:
             disp_dataset = self._phonon_dataset
 
+        if not forces_in_dataset(disp_dataset):
+            raise RuntimeError("Forces are not set in the dataset.")
+
         if is_compact_fc:
             p2s_map = self._phonon_primitive.p2s_map
         else:
@@ -1548,10 +1554,7 @@ def produce_fc2(
             )
         else:
             if "displacements" in disp_dataset:
-                msg = (
-                    "fc_calculator has to be set to produce force "
-                    "constans from this dataset for fc2."
-                )
+                msg = "fc_calculator to solve fc2 has to be set."
                 raise RuntimeError(msg)
             self._fc2 = get_phonopy_fc2(
                 self._phonon_supercell,

phono3py/cui/create_force_constants.py

@@ -217,51 +217,40 @@ def parse_forces(
     dataset = None
     calculator = phono3py.calculator
 
+    if phono3py.phonon_supercell is None or fc_type == "fc3":
+        natom = len(phono3py.supercell)
+    else:
+        natom = len(phono3py.phonon_supercell)
+
     # Get dataset from ph3py_yaml. dataset can be None.
     # physical_units can be overwritten if calculator is found in ph3py_yaml.
     if ph3py_yaml:
         dataset = _extract_dataset_from_ph3py_yaml(ph3py_yaml, fc_type)
-    if dataset and ph3py_yaml.calculator:
-        calculator = ph3py_yaml.calculator
+        if dataset:
+            filename_read_from = phono3py_yaml_filename
+            if ph3py_yaml.calculator:
+                calculator = ph3py_yaml.calculator
 
     physical_units = get_default_physical_units(calculator)
 
-    if phono3py.phonon_supercell is None or fc_type == "fc3":
-        natom = len(phono3py.supercell)
-    else:
-        natom = len(phono3py.phonon_supercell)
-
-    if dataset:
-        filename_read_from = phono3py_yaml_filename
-
-        # Units of displacements and forces are converted. If forces don't
-        # exist, the convesion will not be performed for forces.
-        if calculator is not None:
-            _convert_unit_in_dataset(
-                dataset,
-                distance_to_A=physical_units["distance_to_A"],
-                force_to_eVperA=physical_units["force_to_eVperA"],
+    # Forces are not yet found in dataset. Then try to read from FORCES_FC3 or
+    # FORCES_FC2.
+    if force_filename is not None:
+        if dataset is None or (dataset is not None and not forces_in_dataset(dataset)):
+            dataset = _read_FORCES_FC3_or_FC2(
+                natom, dataset, fc_type, filename=force_filename, log_level=log_level
             )
-
-    # Try to read FORCES_FC* if type-2 and return dataset.
-    # None is returned unless type-2.
-    # can emit FileNotFoundError.
-    if dataset is None or (dataset is not None and not forces_in_dataset(dataset)):
-        _dataset = read_type2_dataset(
-            natom, filename=force_filename, log_level=log_level
+            if dataset:
+                filename_read_from = force_filename
+
+    # Units of displacements and forces are converted. If forces don't
+    # exist, the convesion will not be performed for forces.
+    if calculator is not None:
+        _convert_unit_in_dataset(
+            dataset,
+            distance_to_A=physical_units["distance_to_A"],
+            force_to_eVperA=physical_units["force_to_eVperA"],
         )
-        # Do not overwrite dataset when _dataset is None.
-        if _dataset:
-            filename_read_from = force_filename
-            dataset = _dataset
-
-            # Units of displacements and forces are converted.
-            if calculator is not None:
-                _convert_unit_in_dataset(
-                    dataset,
-                    distance_to_A=physical_units["distance_to_A"],
-                    force_to_eVperA=physical_units["force_to_eVperA"],
-                )
 
     assert dataset is not None
 
@@ -286,29 +275,6 @@ def parse_forces(
             if log_level:
                 print("Cutoff-pair-distance: %f" % cutoff_pair_distance)
 
-    # Type-1 FORCES_FC*.
-    # dataset comes either from disp_fc*.yaml or phono3py*.yaml.
-    if not forces_in_dataset(dataset):
-        if force_filename is not None:
-            if fc_type == "fc3":
-                parse_FORCES_FC3(dataset, filename=force_filename)
-            else:
-                parse_FORCES_FC2(dataset, filename=force_filename)
-
-            if log_level:
-                print(
-                    f'Sets of supercell forces were read from "{force_filename}".',
-                    flush=True,
-                )
-
-        # Unit of displacements is already converted.
-        # Therefore, only unit of forces is converted.
-        if calculator is not None:
-            _convert_unit_in_dataset(
-                dataset,
-                force_to_eVperA=physical_units["force_to_eVperA"],
-            )
-
     return dataset
 
 
@@ -437,20 +403,40 @@ def _read_phono3py_fc2(phono3py, symmetrize_fc2, input_filename, log_level):
     phono3py.fc2 = phonon_fc2
 
 
-def read_type2_dataset(natom, filename="FORCES_FC3", log_level=0) -> Optional[dict]:
-    """Read type-2 FORCES_FC3."""
+def _read_FORCES_FC3_or_FC2(
+    natom: int,
+    dataset: Optional[dict],
+    fc_type: str,
+    filename: str = "FORCES_FC3",
+    log_level: int = 0,
+) -> Optional[dict]:
+    """Read FORCES_FC3 or FORCES_FC2.
+
+    Read the first line of forces file to determine the type of the file.
+
+    """
     if filename is None or not pathlib.Path(filename).exists():
         return None
 
     with open(filename, "r") as f:
         len_first_line = get_length_of_first_line(f)
-        if len_first_line == 6:
-            dataset = get_dataset_type2(f, natom)
+        if len_first_line == 6:  # Type-2
+            _dataset = get_dataset_type2(f, natom)
             if log_level:
-                n_disp = len(dataset["displacements"])
+                n_disp = len(_dataset["displacements"])
                 print(f'{n_disp} snapshots were found in "{filename}".')
-        else:
-            dataset = None
+            return _dataset
+
+    # Type-1
+    if fc_type == "fc3":
+        parse_FORCES_FC3(dataset, filename)
+    else:
+        parse_FORCES_FC2(dataset, filename)
+    if log_level:
+        print(
+            f'Sets of supercell forces were read from "{filename}".',
+            flush=True,
+        )
     return dataset