diff --git a/src/nomad_simulations/schema_packages/outputs.py b/src/nomad_simulations/schema_packages/outputs.py
index abc8570be..d6e20c869 100644
--- a/src/nomad_simulations/schema_packages/outputs.py
+++ b/src/nomad_simulations/schema_packages/outputs.py
@@ -1,7 +1,6 @@
 from typing import TYPE_CHECKING, Optional
 
 import numpy as np
-from nomad.datamodel.data import ArchiveSection
 from nomad.metainfo import Quantity, SubSection
 
 if TYPE_CHECKING:
@@ -22,7 +21,6 @@
     ElectronicEigenvalues,
     ElectronicGreensFunction,
     ElectronicSelfEnergy,
-    FermiLevel,
     FermiSurface,
     HoppingMatrix,
     HybridizationFunction,
@@ -68,8 +66,6 @@ class Outputs(Time):
         """,
     )
 
-    fermi_levels = SubSection(sub_section=FermiLevel.m_def, repeats=True)
-
     chemical_potentials = SubSection(sub_section=ChemicalPotential.m_def, repeats=True)
 
     crystal_field_splittings = SubSection(
diff --git a/src/nomad_simulations/schema_packages/physical_property.py b/src/nomad_simulations/schema_packages/physical_property.py
index 3d56820f8..aa88f3bad 100644
--- a/src/nomad_simulations/schema_packages/physical_property.py
+++ b/src/nomad_simulations/schema_packages/physical_property.py
@@ -1,73 +1,30 @@
-from functools import wraps
-from typing import TYPE_CHECKING, Optional
+from typing import TYPE_CHECKING
 
-import numpy as np
 from nomad import utils
-from nomad.datamodel.data import ArchiveSection
 from nomad.datamodel.metainfo.basesections.v2 import Entity
-from nomad.metainfo import URL, MEnum, Quantity, Reference, SectionProxy, SubSection
+from nomad.datamodel.metainfo.plot import PlotlyFigure, PlotSection
+from nomad.metainfo import URL, Quantity, Reference, SectionProxy, SubSection
 
-if TYPE_CHECKING:
-    from nomad.datamodel.datamodel import EntryArchive
-    from nomad.metainfo import Context, Section
-    from structlog.stdlib import BoundLogger
-
-from nomad_simulations.schema_packages.model_method import BaseModelMethod
 from nomad_simulations.schema_packages.numerical_settings import SelfConsistency
-from nomad_simulations.schema_packages.variables import Variables
 
-# We add `logger` for the `validate_quantity_wrt_value` decorator
 logger = utils.get_logger(__name__)
 
-
-def validate_quantity_wrt_value(name: str = ''):
-    """
-    Decorator to validate the existence of a quantity and its shape with respect to the `PhysicalProperty.value`
-    before calling a method. An example can be found in the module `properties/band_structure.py` for the method
-    `ElectronicEigenvalues.order_eigenvalues()`.
-
-    Args:
-        name (str, optional): The name of the `quantity` to validate. Defaults to ''.
-    """
-
-    def decorator(func):
-        @wraps(func)
-        def wrapper(self, *args, **kwargs):
-            # Checks if `quantity` is defined
-            quantity = getattr(self, name, None)
-            if quantity is None or len(quantity) == 0:
-                logger.warning(f'The quantity `{name}` is not defined.')
-                return False
-
-            # Checks if `value` exists and has the same shape as `quantity`
-            value = getattr(self, 'value', None)
-            if value is None:
-                logger.warning('The quantity `value` is not defined.')
-                return False
-            if value is not None and value.shape != quantity.shape:
-                logger.warning(
-                    f'The shape of the quantity `{name}` does not match the shape of the `value`.'
-                )
-                return False
-
-            return func(self, *args, **kwargs)
-
-        return wrapper
-
-    return decorator
+if TYPE_CHECKING:
+    from nomad.datamodel.metainfo import BoundLogger
 
 
-class PhysicalProperty(ArchiveSection):
+class PhysicalProperty(PlotSection):
     """
-    A base section used to define the physical properties obtained in a simulation, experiment, or in a post-processing
-    analysis. The main quantity of the `PhysicalProperty` is `value`, whose instantiation has to be overwritten in the derived classes
-    when inheriting from `PhysicalProperty`. It contains `variables`, to define the variables over which the physical property varies (see variables.py).
-    This class can also store several string identifiers and quantities for referencing and establishing the character of a physical property.
+    A base section for computational output properties, containing all relevant
+    (meta)data. This includes support for visualization and plotting.
+
+    - Supports the definition and use of `value` for the main property data.
+    - Allows for the inclusion of contributions (e.g., via the `contribution_type` attribute and
+      possible subsections), enabling representation of properties that are composed of multiple
+      parts or sources.
+    - Inherits from `PlotSection`, enabling direct integration with plotting and visualization tools.
     """
 
-    # TODO add `errors`
-    # TODO add `smearing`
-
     name = Quantity(
         type=str,
         description="""
@@ -77,28 +34,28 @@ class PhysicalProperty(ArchiveSection):
 
     iri = Quantity(
         type=URL,
+        default='',
         description="""
-        Internationalized Resource Identifier (IRI) of the physical property defined in the FAIRmat
-        taxonomy, https://fairmat-nfdi.github.io/fairmat-taxonomy/.
+        Internationalized Resource Identifier (IRI) pointing to a definition,
+        typically within a larger, ontological framework.
         """,
     )
 
-    source = Quantity(
-        type=MEnum('simulation', 'measurement', 'analysis'),
-        default='simulation',
+    type = Quantity(
+        type=str,
         description="""
-        Source of the physical property. This quantity is related with the `Activity` performed to obtain the physical
-        property. Example: an `ElectronicBandGap` can be obtained from a `'simulation'` or in a `'measurement'`.
+        Type categorization of the physical property. Example: an `ElectronicBandGap` can be `'direct'`
+        or `'indirect'`.
         """,
     )
 
-    type = Quantity(
+    contribution_type = Quantity(
         type=str,
+        default=None,
         description="""
-        Type categorization of the physical property. Example: an `ElectronicBandGap` can be `'direct'`
-        or `'indirect'`.
+        Type of contribution to the physical property. Hence, only applies to `contributions` instances.
+        Example: `TotalEnergy` may have contributions like _kinetic_, _potential_, etc.
         """,
-        # ! add more examples in the description to improve the understanding of this quantity
     )
 
     label = Quantity(
@@ -107,13 +64,8 @@ class PhysicalProperty(ArchiveSection):
         Label for additional classification of the physical property. Example: an `ElectronicBandGap`
         can be labeled as `'DFT'` or `'GW'` depending on the methodology used to calculate it.
         """,
-        # ! add more examples in the description to improve the understanding of this quantity
-    )
-
-    # variables = SubSection(sub_section=Variables.m_def, repeats=True)
+    )  # TODO: specify use better
 
-    # * `value` must be overwritten in the derived classes defining its type, unit, and description
-    # TODO use abstract to enforce policy?
     value: Quantity = None
 
     entity_ref = Quantity(
@@ -124,15 +76,7 @@ class PhysicalProperty(ArchiveSection):
             cell. In the first case, `outputs.model_system_ref` (see outputs.py) will point to the `ModelSystem` section,
             while in the second case, `entity_ref` will point to `AtomsState` section (see atoms_state.py).
         """,
-    )
-
-    physical_property_ref = Quantity(
-        type=Reference(SectionProxy('PhysicalProperty')),
-        description="""
-        Reference to the `PhysicalProperty` section from which the physical property was derived. If `physical_property_ref`
-        is populated, the quantity `is_derived` is set to True via normalization.
-        """,
-    )
+    )  # TODO: only used for electronic states, remove
 
     is_derived = Quantity(
         type=bool,
@@ -146,13 +90,21 @@ class PhysicalProperty(ArchiveSection):
         """,
     )
 
+    physical_property_ref = Quantity(
+        type=Reference(SectionProxy('PhysicalProperty')),
+        description="""
+        Reference to the `PhysicalProperty` section from which the physical property was derived. If `physical_property_ref`
+        is populated, the quantity `is_derived` is set to True via normalization.
+        """,
+    )
+
     is_scf_converged = Quantity(
         type=bool,
         description="""
         Flag indicating whether the physical property is converged or not after a SCF process. This quantity is connected
         with `SelfConsistency` defined in the `numerical_settings.py` module.
         """,
-    )
+    )  # ? tie to calculation, not individual property
 
     self_consistency_ref = Quantity(
         type=SelfConsistency,
@@ -160,95 +112,170 @@ class PhysicalProperty(ArchiveSection):
         Reference to the `SelfConsistency` section that defines the numerical settings to converge the
         physical property (see numerical_settings.py).
         """,
+    )  # ? remove
+
+    contributions = SubSection(
+        section_def=SectionProxy('PhysicalProperty'),
+        repeats=True,
+        description="""
+        Shallow list of contributions to the physical property.
+        Does not necessarily entail a (full) partioning.
+        """,
     )
+    # TODO: would be wishful to have `section_def` be a stripped down version of PhysicalProperty
+    # that gets automatically updated when extending PhysicalProperty
+    # should be discussed with @TLCFEM
+
+    def _is_derived(self) -> bool:
+        """
+        Resolves whether the physical property is derived or not.
+
+        Returns:
+            (bool): The flag indicating whether the physical property is derived or not.
+        """
+        return self.physical_property_ref is not None
+
+    def _is_contribution(self) -> bool:
+        """
+        Determines if this instance is a contribution by checking if it's contained
+        in a parent's contributions list.
+
+        Returns:
+            (bool): True if this instance is a contribution, False otherwise.
+        """
+        if hasattr(self, 'm_parent') and self.m_parent:
+            parent_section = self.m_parent
+            # If parent has contributions containing this instance, we are a contribution
+            if (
+                hasattr(parent_section, 'contributions')
+                and parent_section.contributions
+            ):
+                if self in parent_section.contributions:
+                    return True
+        return False
+
+    def _validate_contributions_structure(self, logger: 'BoundLogger') -> bool:
+        """
+        Validates that contributions do not contain nested contributions.
+        This prevents recursive contribution structures which are not intended.
+        Only runs for top-level PhysicalProperty instances, not for contributions themselves.
+
+        Args:
+            logger: Logger instance for error reporting.
 
-    # @property
-    # def variables_shape(self) -> Optional[list]:
-    #    """
-    #    Shape of the variables over which the physical property varies. This is extracted from
-    #    `Variables.n_points` and appended in a list.
-    #
-    #    Example, a physical property which varies with `Temperature` and `ElectricField` will
-    #    return `variables_shape = [n_temperatures, n_electric_fields]`.
-    #
-    #    Returns:
-    #        (list): The shape of the variables over which the physical property varies.
-    #    """
-    #    if self.variables is not None:
-    #        return [v.get_n_points(logger) for v in self.variables]
-    #    return []
-
-    @property
-    def full_shape(self) -> list[int]:
+        Returns:
+            (bool): True if validation passes, False if nested contributions are found.
         """
-        Full shape of the physical property. This quantity is calculated as a concatenation of the `variables_shape`
-        and `value.shape`:
+        # Skip validation for contribution instances
+        if self._is_contribution():
+            return True
+
+        if not self.contributions:
+            return True
+
+        has_nested_contributions = False
+        for i, contribution in enumerate(self.contributions):
+            if hasattr(contribution, 'contributions') and contribution.contributions:
+                logger.error(
+                    f'Contribution {i} in {self.__class__.__name__} contains nested contributions. '
+                    'Contributions should not have their own contributions subsection populated.'
+                )
+                has_nested_contributions = True
 
-            `full_shape = variables_shape + value.shape`
+        return not has_nested_contributions
 
-        Example: a physical property which is a 3D vector and varies with `variables=[Temperature, ElectricField]`
-        will have `value.shape=[3]`, `variables_shape=[n_temperatures, n_electric_fields]`, and thus
-        `full_shape=[n_temperatures, n_electric_fields, 3]`.
+    def _validate_contribution_type(self, logger) -> bool:
+        """
+        Validates that contribution_type is only set for contribution instances
+        and is not set for top-level PhysicalProperty instances.
+
+        Args:
+            logger: Logger instance for error reporting.
 
         Returns:
-            (list): The full shape of the physical property.
+            (bool): True if validation passes, False if contribution_type is incorrectly set.
         """
-        value_shape = self.value.shape if self.value is not None else []
-        # return self.variables_shape + value_shape
-        return [] + value_shape
-
-    def __init__(
-        self, m_def: 'Section' = None, m_context: 'Context' = None, **kwargs
-    ) -> None:
-        super().__init__(m_def, m_context, **kwargs)
-        # Checking if IRI is defined
-        if self.iri is None:
-            logger.warning(
-                'The used property is not defined in the FAIRmat taxonomy (https://fairmat-nfdi.github.io/fairmat-taxonomy/). You can contribute there if you want to extend the list of available materials properties.'
+        is_contribution = self._is_contribution()
+
+        # Check for incorrect usage
+        if not is_contribution and self.contribution_type is not None:
+            logger.error(
+                f'{self.__class__.__name__} has contribution_type set but is not a contribution. '
+                'contribution_type should only be set for instances in the contributions subsection.'
             )
+            return False
 
-    def __setattr__(self, name, value):
-        if name == 'value':
-            try:
-                value = np.array(value)
-                # self.__class__.value.shape = ['*'] * value.ndim
-            except Exception:
-                pass
-        return super().__setattr__(name, value)
+        return True
 
-    def _is_derived(self) -> bool:
+    def plot(self, **kwargs) -> list[PlotlyFigure]:
         """
-        Resolves whether the physical property is derived or not.
+        Placeholder for a method to plot the physical property. This method should be overridden in derived classes
+        to provide specific plotting functionality.
 
         Returns:
-            (bool): The flag indicating whether the physical property is derived or not.
+            (list[PlotlyFigure]): A list of PlotlyFigure objects representing the physical property.
         """
-        if self.physical_property_ref is not None:
-            return True
-        return False
+        return []
 
-    def normalize(self, archive: 'EntryArchive', logger: 'BoundLogger') -> None:
-        self.is_derived = self._is_derived()
+    def sub_plots(self, **kwargs) -> None:
+        """
+        Collects plots from `self.contributions` and overlays them onto the target figure.
+        """
+        if not self.contributions or not self.figures:
+            return
+
+        try:
+            target_figure = self.figures[kwargs.get('target_indices', -1)]
+        except (IndexError, TypeError):
+            return
+
+        if target_figure.figure:
+            figure_dict = target_figure.figure.copy()
+        else:
+            figure_dict = {'data': [], 'layout': {}}
+
+        for contribution in self.contributions:
+            # Use existing figures if already normalized, otherwise call plot()
+            plots = (
+                contribution.figures
+                if contribution.figures
+                else contribution.plot(**kwargs)
+            )
 
+            if plots:
+                for plot in plots:
+                    if hasattr(plot, 'figure') and plot.figure:
+                        plot_data = plot.figure.get('data', [])
+                        for trace in plot_data:
+                            figure_dict['data'].append(trace)
 
-class PropertyContribution(PhysicalProperty):
-    """
-    Abstract physical property section linking a property contribution to a contribution
-    from some method.
+        target_figure.figure = figure_dict
 
-    Abstract class for incorporating specific contributions of a physical property, while
-    linking this contribution to a specific component (of class `BaseModelMethod`) of the
-    over `ModelMethod` using the `model_method_ref` quantity.
-    """
+    def normalize(self, *args, **kwargs) -> None:
+        # check whether already normalized
+        if self.m_cache.get('_is_normalized', False):
+            return
+        else:
+            self.m_cache['_is_normalized'] = True
 
-    model_method_ref = Quantity(
-        type=BaseModelMethod,
-        description="""
-        Reference to the `ModelMethod` section to which the property is linked to.
-        """,
-    )
+        # perform own normalization
+        super().normalize(*args, **kwargs)
+
+        self.is_derived = self._is_derived()
+
+        # validate contributions structure and contribution_type usage
+        logger_arg = args[1] if len(args) > 1 else logger
+        self._validate_contributions_structure(logger_arg)
+        self._validate_contribution_type(logger_arg)
+
+        for contribution in self.contributions:
+            if hasattr(contribution, 'normalize'):
+                contribution.normalize(*args, **kwargs)
+
+        if plot_figures := self.plot(**kwargs):
+            self.figures.extend(plot_figures)
+        self.sub_plots(**kwargs)
 
-    def normalize(self, archive, logger) -> None:
-        super().normalize(archive, logger)
-        if not self.name:
-            self.name = self.get('model_method_ref', {}).get('name')
+        # set names last, they may depend other normalized properties
+        if self.m_def.name is not None:
+            self.name = self.m_def.name
diff --git a/src/nomad_simulations/schema_packages/properties/__init__.py b/src/nomad_simulations/schema_packages/properties/__init__.py
index 69470a4ff..de70abb90 100644
--- a/src/nomad_simulations/schema_packages/properties/__init__.py
+++ b/src/nomad_simulations/schema_packages/properties/__init__.py
@@ -1,14 +1,12 @@
 from .band_gap import ElectronicBandGap
 from .band_structure import ElectronicBandStructure, ElectronicEigenvalues, Occupancy
 from .energies import (
-    EnergyContribution,
-    FermiLevel,
     KineticEnergy,
     PotentialEnergy,
     TotalEnergy,
 )
 from .fermi_surface import FermiSurface
-from .forces import BaseForce, ForceContribution, TotalForce
+from .forces import TotalForce
 from .greens_function import (
     ElectronicGreensFunction,
     ElectronicSelfEnergy,
diff --git a/src/nomad_simulations/schema_packages/properties/band_gap.py b/src/nomad_simulations/schema_packages/properties/band_gap.py
index 5e97acf0c..2667aa107 100644
--- a/src/nomad_simulations/schema_packages/properties/band_gap.py
+++ b/src/nomad_simulations/schema_packages/properties/band_gap.py
@@ -57,12 +57,6 @@ class ElectronicBandGap(PhysicalProperty):
         """,
     )
 
-    def __init__(
-        self, m_def: 'Section' = None, m_context: 'Context' = None, **kwargs
-    ) -> None:
-        super().__init__(m_def, m_context, **kwargs)
-        self.name = self.m_def.name
-
     def resolve_type(self, logger: 'BoundLogger') -> Optional[str]:
         """
         Resolves the `type` of the electronic band gap based on the stored `momentum_transfer` values.
diff --git a/src/nomad_simulations/schema_packages/properties/band_structure.py b/src/nomad_simulations/schema_packages/properties/band_structure.py
index 30c81955c..f03c9cd4a 100644
--- a/src/nomad_simulations/schema_packages/properties/band_structure.py
+++ b/src/nomad_simulations/schema_packages/properties/band_structure.py
@@ -9,15 +9,11 @@
 
 if TYPE_CHECKING:
     from nomad.datamodel.datamodel import EntryArchive
-    from nomad.metainfo import Context, Section
     from structlog.stdlib import BoundLogger
 
 from nomad_simulations.schema_packages.atoms_state import AtomsState, OrbitalsState
 from nomad_simulations.schema_packages.numerical_settings import KSpace
-from nomad_simulations.schema_packages.physical_property import (
-    PhysicalProperty,
-    validate_quantity_wrt_value,
-)
+from nomad_simulations.schema_packages.physical_property import PhysicalProperty
 from nomad_simulations.schema_packages.properties.band_gap import ElectronicBandGap
 from nomad_simulations.schema_packages.properties.fermi_surface import FermiSurface
 from nomad_simulations.schema_packages.utils import get_sibling_section
@@ -32,8 +28,6 @@ class BaseElectronicEigenvalues(PhysicalProperty):
     A base section used to define basic quantities for the `ElectronicEigenvalues`  and `ElectronicBandStructure` properties.
     """
 
-    iri = ''
-
     n_bands = Quantity(
         type=np.int32,
         description="""
@@ -118,22 +112,21 @@ class ElectronicEigenvalues(BaseElectronicEigenvalues):
         """,
     )
 
-    def __init__(
-        self, m_def: 'Section' = None, m_context: 'Context' = None, **kwargs
-    ) -> None:
-        super().__init__(m_def, m_context, **kwargs)
-        self.name = self.m_def.name
-
-    @validate_quantity_wrt_value(name='occupation')
-    def order_eigenvalues(self) -> Union[bool, tuple[pint.Quantity, np.ndarray]]:
+    def order_eigenvalues(self) -> tuple[pint.Quantity, np.ndarray] | tuple[()]:
         """
         Order the eigenvalues based on the `value` and `occupation`. The return `value` and
         `occupation` are flattened.
 
         Returns:
-            (Union[bool, tuple[pint.Quantity, np.ndarray]]): The flattened and sorted `value` and `occupation`. If validation
-            fails, then it returns `False`.
+            (tuple[pint.Quantity, np.ndarray] | tuple[()]): The flattened and sorted `value` and `occupation`. If validation
+            fails, then it returns an empty tuple.
         """
+        # Validation: check if both value and occupation exist and have same shape
+        if self.value is None or self.occupation is None:
+            return ()
+        if self.value.shape != self.occupation.shape:
+            return ()
+
         total_shape = np.prod(self.value.shape)
 
         # Order the indices in the flattened list of `value`
@@ -164,13 +157,14 @@ def resolve_homo_lumo_eigenvalues(
             `lowest_unoccupied` eigenvalues.
         """
         # Sorting `value` and `occupation`
-        if not self.order_eigenvalues():  # validation fails
+        if ordered_results := self.order_eigenvalues():
+            sorted_value, sorted_occupation = ordered_results
+            sorted_value_unit = sorted_value.u
+            sorted_value = sorted_value.magnitude
+        else:
             if self.highest_occupied is not None and self.lowest_unoccupied is not None:
                 return self.highest_occupied, self.lowest_unoccupied
             return None, None
-        sorted_value, sorted_occupation = self.order_eigenvalues()
-        sorted_value_unit = sorted_value.u
-        sorted_value = sorted_value.magnitude
 
         # Binary search ot find the transition point between `occupation = 2` and `occupation = 0`
         homo = self.highest_occupied
@@ -309,12 +303,6 @@ class ElectronicBandStructure(ElectronicEigenvalues):
 
     k_path = SubSection(sub_section=KLinePath.m_def)
 
-    def __init__(
-        self, m_def: 'Section' = None, m_context: 'Context' = None, **kwargs
-    ) -> None:
-        super().__init__(m_def, m_context, **kwargs)
-        self.name = self.m_def.name
-
 
 class Occupancy(PhysicalProperty):
     """
@@ -357,10 +345,4 @@ class Occupancy(PhysicalProperty):
         """,
     )
 
-    def __init__(
-        self, m_def: 'Section' = None, m_context: 'Context' = None, **kwargs
-    ) -> None:
-        super().__init__(m_def, m_context, **kwargs)
-        self.name = self.m_def.name
-
     # TODO add extraction from `ElectronicEigenvalues.occupation`
diff --git a/src/nomad_simulations/schema_packages/properties/energies.py b/src/nomad_simulations/schema_packages/properties/energies.py
index 833196b28..dc57f5b59 100644
--- a/src/nomad_simulations/schema_packages/properties/energies.py
+++ b/src/nomad_simulations/schema_packages/properties/energies.py
@@ -1,17 +1,7 @@
-from typing import TYPE_CHECKING
-
 import numpy as np
-from nomad.metainfo import Context, Quantity, Section, SubSection
-
-if TYPE_CHECKING:
-    from nomad.datamodel.datamodel import EntryArchive
-    from nomad.metainfo import Context, Section
-    from structlog.stdlib import BoundLogger
+from nomad.metainfo import Quantity, SectionProxy, SubSection
 
-from nomad_simulations.schema_packages.physical_property import (
-    PhysicalProperty,
-    PropertyContribution,
-)
+from nomad_simulations.schema_packages.physical_property import PhysicalProperty
 
 ##################
 # Abstract classes
@@ -33,42 +23,11 @@ class BaseEnergy(PhysicalProperty):
     )
 
 
-class EnergyContribution(BaseEnergy, PropertyContribution):
-    """
-    Abstract class for incorporating specific energy contributions to the `TotalEnergy`.
-    The inheritance from `PropertyContribution` allows to link this contribution to a
-    specific component (of class `BaseModelMethod`) of the over `ModelMethod` using the
-    `model_method_ref` quantity.
-
-    For example, for a force field calculation, the `model_method_ref` may point to a
-    particular potential type (e.g., a Lennard-Jones potential between atom types X and Y),
-    while for a DFT calculation, it may point to a particular electronic interaction term
-    (e.g., 'XC' for the exchange-correlation term, or 'Hartree' for the Hartree term).
-    Then, the contribution will be named according to this model component and the `value`
-    quantity will contain the energy contribution from this component evaluated over all
-    relevant atoms or electrons or as a function of them.
-    """
-
-
 ####################################
 # List of specific energy properties
 ####################################
 
 
-class FermiLevel(BaseEnergy):
-    """
-    Energy required to add or extract a charge from a material at zero temperature. It can be also defined as the chemical potential at zero temperature.
-    """
-
-    iri = 'http://fairmat-nfdi.eu/taxonomy/FermiLevel'
-
-    def __init__(
-        self, m_def: 'Section' = None, m_context: 'Context' = None, **kwargs
-    ) -> None:
-        super().__init__(m_def, m_context, **kwargs)
-        self.name = self.m_def.name
-
-
 #! The only issue with this structure is that total energy will never be a sum of its contributions,
 #! since kinetic energy lives separately, but I think maybe this is ok?
 class TotalEnergy(BaseEnergy):
@@ -78,13 +37,7 @@ class TotalEnergy(BaseEnergy):
     """
 
     # ? add a generic contributions quantity to PhysicalProperty
-    contributions = SubSection(sub_section=EnergyContribution.m_def, repeats=True)
-
-    def __init__(
-        self, m_def: 'Section' = None, m_context: 'Context' = None, **kwargs
-    ) -> None:
-        super().__init__(m_def, m_context, **kwargs)
-        self.name = self.m_def.name
+    contributions = SubSection(sub_section=SectionProxy('BaseEnergy'), repeats=True)
 
 
 # ? Separate quantities for nuclear and electronic KEs?
@@ -93,20 +46,8 @@ class KineticEnergy(BaseEnergy):
     Physical property section describing the kinetic energy of a (sub)system.
     """
 
-    def __init__(
-        self, m_def: 'Section' = None, m_context: 'Context' = None, **kwargs
-    ) -> None:
-        super().__init__(m_def, m_context, **kwargs)
-        self.name = self.m_def.name
-
 
 class PotentialEnergy(BaseEnergy):
     """
     Physical property section describing the potential energy of a (sub)system.
     """
-
-    def __init__(
-        self, m_def: 'Section' = None, m_context: 'Context' = None, **kwargs
-    ) -> None:
-        super().__init__(m_def, m_context, **kwargs)
-        self.name = self.m_def.name
diff --git a/src/nomad_simulations/schema_packages/properties/fermi_surface.py b/src/nomad_simulations/schema_packages/properties/fermi_surface.py
index c1defbfa1..51dd28fc6 100644
--- a/src/nomad_simulations/schema_packages/properties/fermi_surface.py
+++ b/src/nomad_simulations/schema_packages/properties/fermi_surface.py
@@ -28,12 +28,5 @@ class FermiSurface(PhysicalProperty):
         """,
     )
 
-    def __init__(
-        self, m_def: 'Section' = None, m_context: 'Context' = None, **kwargs
-    ) -> None:
-        super().__init__(m_def, m_context, **kwargs)
-        # ! `n_bands` need to be set up during initialization of the class
-        self.name = self.m_def.name
-
     def normalize(self, archive: 'EntryArchive', logger: 'BoundLogger') -> None:
         super().normalize(archive, logger)
diff --git a/src/nomad_simulations/schema_packages/properties/forces.py b/src/nomad_simulations/schema_packages/properties/forces.py
index 5d4a5bc33..c3d9aab26 100644
--- a/src/nomad_simulations/schema_packages/properties/forces.py
+++ b/src/nomad_simulations/schema_packages/properties/forces.py
@@ -1,24 +1,14 @@
-from typing import TYPE_CHECKING
-
 import numpy as np
-from nomad.metainfo import Context, Quantity, Section, SubSection
-
-if TYPE_CHECKING:
-    from nomad.datamodel.datamodel import EntryArchive
-    from nomad.metainfo import Context, Section
-    from structlog.stdlib import BoundLogger
+from nomad.metainfo import Quantity, SectionProxy, SubSection
 
-from nomad_simulations.schema_packages.physical_property import (
-    PhysicalProperty,
-    PropertyContribution,
-)
+from nomad_simulations.schema_packages.physical_property import PhysicalProperty
 
 ##################
 # Abstract classes
 ##################
 
 
-class BaseForce(PhysicalProperty):
+class TotalForce(PhysicalProperty):
     """
     Abstract class used to define a common `value` quantity with the appropriate units
     for different types of forces, which avoids repeating the definitions for each
@@ -33,39 +23,11 @@ class BaseForce(PhysicalProperty):
         """,
     )
 
-
-class ForceContribution(BaseForce, PropertyContribution):
-    """
-    Abstract class for incorporating specific force contributions to the `TotalForce`.
-    The inheritance from `PropertyContribution` allows to link this contribution to a
-    specific component (of class `BaseModelMethod`) of the over `ModelMethod` using the
-    `model_method_ref` quantity.
-
-    For example, for a force field calculation, the `model_method_ref` may point to a
-    particular potential type (e.g., a Lennard-Jones potential between atom types X and Y),
-    while for a DFT calculation, it may point to a particular electronic interaction term
-    (e.g., 'XC' for the exchange-correlation term, or 'Hartree' for the Hartree term).
-    Then, the contribution will be named according to this model component and the `value`
-    quantity will contain the force contribution from this component evaluated over all
-    relevant atoms or electrons or as a function of them.
-    """
-
-
-###################################
-# List of specific force properties
-###################################
-
-
-class TotalForce(BaseForce):
-    """
-    The total force on a system. `contributions` specify individual force
-    contributions to the `TotalForce`.
-    """
-
-    contributions = SubSection(sub_section=ForceContribution.m_def, repeats=True)
-
-    def __init__(
-        self, m_def: 'Section' = None, m_context: 'Context' = None, **kwargs
-    ) -> None:
-        super().__init__(m_def, m_context, **kwargs)
-        self.name = self.m_def.name
+    contributions = SubSection(
+        sub_section=SectionProxy('TotalForce'),
+        repeats=True,
+        description="""
+        The contributions to the total force. Each contribution is a specific force
+        component, such as the force from a specific potential or interaction.
+        """,
+    )
diff --git a/src/nomad_simulations/schema_packages/properties/greens_function.py b/src/nomad_simulations/schema_packages/properties/greens_function.py
index bca9cda62..c9fae7828 100644
--- a/src/nomad_simulations/schema_packages/properties/greens_function.py
+++ b/src/nomad_simulations/schema_packages/properties/greens_function.py
@@ -189,12 +189,6 @@ class ElectronicGreensFunction(BaseGreensFunction):
         """,
     )
 
-    def __init__(
-        self, m_def: 'Section' = None, m_context: 'Context' = None, **kwargs
-    ) -> None:
-        super().__init__(m_def, m_context, **kwargs)
-        self.name = self.m_def.name
-
 
 class ElectronicSelfEnergy(BaseGreensFunction):
     """
@@ -211,12 +205,6 @@ class ElectronicSelfEnergy(BaseGreensFunction):
         """,
     )
 
-    def __init__(
-        self, m_def: 'Section' = None, m_context: 'Context' = None, **kwargs
-    ) -> None:
-        super().__init__(m_def, m_context, **kwargs)
-        self.name = self.m_def.name
-
 
 class HybridizationFunction(BaseGreensFunction):
     """
@@ -233,12 +221,6 @@ class HybridizationFunction(BaseGreensFunction):
         """,
     )
 
-    def __init__(
-        self, m_def: 'Section' = None, m_context: 'Context' = None, **kwargs
-    ) -> None:
-        super().__init__(m_def, m_context, **kwargs)
-        self.name = self.m_def.name
-
 
 class QuasiparticleWeight(PhysicalProperty):
     """
@@ -320,12 +302,6 @@ class QuasiparticleWeight(PhysicalProperty):
         """,
     )
 
-    def __init__(
-        self, m_def: 'Section' = None, m_context: 'Context' = None, **kwargs
-    ) -> None:
-        super().__init__(m_def, m_context, **kwargs)
-        self.name = self.m_def.name
-
     def resolve_system_correlation_strengths(self) -> str:
         """
         Resolves the `system_correlation_strengths` of the quasiparticle weight based on the stored `value` values.
diff --git a/src/nomad_simulations/schema_packages/properties/hopping_matrix.py b/src/nomad_simulations/schema_packages/properties/hopping_matrix.py
index b5687f29f..445fa4f8f 100644
--- a/src/nomad_simulations/schema_packages/properties/hopping_matrix.py
+++ b/src/nomad_simulations/schema_packages/properties/hopping_matrix.py
@@ -43,12 +43,6 @@ class HoppingMatrix(PhysicalProperty):
         """,
     )
 
-    def __init__(
-        self, m_def: 'Section' = None, m_context: 'Context' = None, **kwargs
-    ) -> None:
-        super().__init__(m_def, m_context, **kwargs)
-        self.name = self.m_def.name
-
     # TODO add normalization to extract DOS, band structure, etc, properties from `HoppingMatrix`
 
 
@@ -75,9 +69,3 @@ class CrystalFieldSplitting(PhysicalProperty):
         at the same Wigner-Seitz point (0, 0, 0).
         """,
     )
-
-    def __init__(
-        self, m_def: 'Section' = None, m_context: 'Context' = None, **kwargs
-    ) -> None:
-        super().__init__(m_def, m_context, **kwargs)
-        self.name = self.m_def.name
diff --git a/src/nomad_simulations/schema_packages/properties/permittivity.py b/src/nomad_simulations/schema_packages/properties/permittivity.py
index 593e29b18..68bb7ce9a 100644
--- a/src/nomad_simulations/schema_packages/properties/permittivity.py
+++ b/src/nomad_simulations/schema_packages/properties/permittivity.py
@@ -26,6 +26,10 @@ class Permittivity(PhysicalProperty):
 
     iri = 'http://fairmat-nfdi.eu/taxonomy/Permittivity'
 
+    # Class-level constants
+    rank = [3, 3]
+    _axes_map = ['xx', 'yy', 'zz']
+
     type = Quantity(
         type=MEnum('static', 'dynamic'),
         description="""
@@ -51,14 +55,6 @@ class Permittivity(PhysicalProperty):
     # ? We need use cases to understand if we need to define contributions to the permittivity tensor.
     # ? `ionic` and `electronic` contributions are common in the literature.
 
-    def __init__(
-        self, m_def: 'Section' = None, m_context: 'Context' = None, **kwargs
-    ) -> None:
-        super().__init__(m_def, m_context, **kwargs)
-        self.rank = [3, 3]
-        self.name = self.m_def.name
-        self._axes_map = ['xx', 'yy', 'zz']
-
     def resolve_type(self) -> str:
         return 'static' if self.frequencies is None else 'dynamic'
 
diff --git a/src/nomad_simulations/schema_packages/properties/spectral_profile.py b/src/nomad_simulations/schema_packages/properties/spectral_profile.py
index 8e0510fdb..abb6384d5 100644
--- a/src/nomad_simulations/schema_packages/properties/spectral_profile.py
+++ b/src/nomad_simulations/schema_packages/properties/spectral_profile.py
@@ -65,11 +65,6 @@ class DOSProfile(SpectralProfile):
         """,
     )
 
-    def __init__(
-        self, m_def: 'Section' = None, m_context: 'Context' = None, **kwargs
-    ) -> None:
-        super().__init__(m_def, m_context, **kwargs)
-
     def resolve_pdos_name(self, logger: 'BoundLogger') -> Optional[str]:
         """
         Resolve the `name` of the projected `DOSProfile` from the `entity_ref` section. This is resolved as:
@@ -191,12 +186,6 @@ class ElectronicDensityOfStates(DOSProfile):
         """,
     )
 
-    def __init__(
-        self, m_def: 'Section' = None, m_context: 'Context' = None, **kwargs
-    ) -> None:
-        super().__init__(m_def, m_context, **kwargs)
-        self.name = self.m_def.name
-
     def resolve_energies_origin(
         self,
         energies_points: pint.Quantity,
@@ -506,13 +495,6 @@ class AbsorptionSpectrum(SpectralProfile):
         """,
     )
 
-    def __init__(
-        self, m_def: 'Section' = None, m_context: 'Context' = None, **kwargs
-    ) -> None:
-        super().__init__(m_def, m_context, **kwargs)
-        # Set the name of the section
-        self.name = self.m_def.name
-
 
 class XASSpectrum(AbsorptionSpectrum):
     """
@@ -535,13 +517,6 @@ class XASSpectrum(AbsorptionSpectrum):
         repeats=False,
     )
 
-    def __init__(
-        self, m_def: 'Section' = None, m_context: 'Context' = None, **kwargs
-    ) -> None:
-        super().__init__(m_def, m_context, **kwargs)
-        # Set the name of the section
-        self.name = self.m_def.name
-
     def generate_from_contributions(self, logger: 'BoundLogger') -> None:
         """
         Generate the `value` of the XAS spectrum by concatenating the XANES and EXAFS contributions. It also concatenates
diff --git a/src/nomad_simulations/schema_packages/properties/thermodynamics.py b/src/nomad_simulations/schema_packages/properties/thermodynamics.py
index 30e3c6c06..88dc020cb 100644
--- a/src/nomad_simulations/schema_packages/properties/thermodynamics.py
+++ b/src/nomad_simulations/schema_packages/properties/thermodynamics.py
@@ -130,15 +130,50 @@ class HelmholtzFreeEnergy(BaseEnergy):
 class ChemicalPotential(BaseEnergy):
     """
     Free energy cost of adding or extracting a particle from a thermodynamic system.
+
+    At finite temperature, the chemical potential determines the equilibrium condition
+    for particle exchange between different phases or subsystems. It can be defined
+    as the partial derivative of the internal energy with respect to particle number
+    at constant entropy and volume, or equivalently as the partial derivative of
+    the Gibbs free energy with respect to particle number at constant temperature
+    and pressure.
     """
 
     iri = 'http://fairmat-nfdi.eu/taxonomy/ChemicalPotential'
 
-    def __init__(
-        self, m_def: 'Section' = None, m_context: 'Context' = None, **kwargs
-    ) -> None:
-        super().__init__(m_def, m_context, **kwargs)
-        self.name = self.m_def.name
+    temperature = Quantity(
+        type=np.float64,
+        unit='kelvin',
+        description="""
+        Temperature at which the chemical potential is calculated.
+        Essential for finite-temperature calculations.
+        """,
+    )
+
+    particle_number = Quantity(
+        type=np.float64,
+        description="""
+        Number of particles (or particle density) for which the chemical potential applies.
+        Can represent electron number, atom number, or other relevant particle count.
+        """,
+    )
+
+    fermi_energy = Quantity(
+        type=np.float64,
+        unit='joule',
+        description="""
+        Fermi energy at T=0K, used as reference for finite-temperature chemical potential.
+        At T=0, the chemical potential equals the Fermi energy.
+        """,
+    )
+
+    type = Quantity(
+        type=str,
+        description="""
+        Type of chemical potential calculation. Examples: 'electronic', 'atomic', 
+        'ionic', 'molecular'. Helps identify what kind of particles this applies to.
+        """,
+    )
 
 
 class HeatCapacity(PhysicalProperty):
@@ -168,12 +203,7 @@ class VirialTensor(BaseEnergy):
     the virial theorem.
     """
 
-    def __init__(
-        self, m_def: 'Section' = None, m_context: 'Context' = None, **kwargs
-    ) -> None:
-        super().__init__(m_def, m_context, **kwargs)
-        self.rank = [3, 3]
-        self.name = self.m_def.name
+    rank = [3, 3]
 
 
 class MassDensity(PhysicalProperty):
@@ -196,16 +226,11 @@ class Hessian(PhysicalProperty):
     describing the local curvature of the energy surface.
     """
 
+    rank = [3, 3]
+
     value = Quantity(
         type=np.float64,
         unit='joule / m ** 2',
         description="""
         """,
     )
-
-    def __init__(
-        self, m_def: 'Section' = None, m_context: 'Context' = None, **kwargs
-    ) -> None:
-        super().__init__(m_def, m_context, **kwargs)
-        self.rank = [3, 3]
-        self.name = self.m_def.name
diff --git a/src/nomad_simulations/schema_packages/variables.py b/src/nomad_simulations/schema_packages/variables.py
index 34e3c8c4d..e485eddfe 100644
--- a/src/nomad_simulations/schema_packages/variables.py
+++ b/src/nomad_simulations/schema_packages/variables.py
@@ -75,6 +75,9 @@ def normalize(self, archive: 'EntryArchive', logger: 'BoundLogger') -> None:
         # Setting `n_points` if these are not defined
         self.n_points = self.get_n_points(logger)
 
+        if self.m_def.name is not None:
+            self.name = self.m_def.name
+
 
 class Temperature(Variables):
     """ """
@@ -88,12 +91,6 @@ class Temperature(Variables):
         """,
     )
 
-    def __init__(
-        self, m_def: 'Section' = None, m_context: 'Context' = None, **kwargs
-    ) -> None:
-        super().__init__(m_def, m_context, **kwargs)
-        self.name = self.m_def.name
-
     def normalize(self, archive: 'EntryArchive', logger: 'BoundLogger') -> None:
         super().normalize(archive, logger)
 
@@ -111,12 +108,6 @@ class Energy2(Variables):
         """,
     )
 
-    def __init__(
-        self, m_def: 'Section' = None, m_context: 'Context' = None, **kwargs
-    ) -> None:
-        super().__init__(m_def, m_context, **kwargs)
-        self.name = self.m_def.name
-
     def normalize(self, archive: 'EntryArchive', logger: 'BoundLogger') -> None:
         super().normalize(archive, logger)
 
@@ -135,12 +126,6 @@ class WignerSeitz(Variables):
         """,
     )
 
-    def __init__(
-        self, m_def: 'Section' = None, m_context: 'Context' = None, **kwargs
-    ) -> None:
-        super().__init__(m_def, m_context, **kwargs)
-        self.name = self.m_def.name
-
     def normalize(self, archive: 'EntryArchive', logger: 'BoundLogger') -> None:
         super().normalize(archive, logger)
 
@@ -157,12 +142,6 @@ class Frequency(Variables):
         """,
     )
 
-    def __init__(
-        self, m_def: 'Section' = None, m_context: 'Context' = None, **kwargs
-    ) -> None:
-        super().__init__(m_def, m_context, **kwargs)
-        self.name = self.m_def.name
-
     def normalize(self, archive: 'EntryArchive', logger: 'BoundLogger') -> None:
         super().normalize(archive, logger)
 
@@ -179,12 +158,6 @@ class MatsubaraFrequency(Variables):
         """,
     )
 
-    def __init__(
-        self, m_def: 'Section' = None, m_context: 'Context' = None, **kwargs
-    ) -> None:
-        super().__init__(m_def, m_context, **kwargs)
-        self.name = self.m_def.name
-
     def normalize(self, archive: 'EntryArchive', logger: 'BoundLogger') -> None:
         super().normalize(archive, logger)
 
@@ -201,12 +174,6 @@ class Time(Variables):
         """,
     )
 
-    def __init__(
-        self, m_def: 'Section' = None, m_context: 'Context' = None, **kwargs
-    ) -> None:
-        super().__init__(m_def, m_context, **kwargs)
-        self.name = self.m_def.name
-
     def normalize(self, archive: 'EntryArchive', logger: 'BoundLogger') -> None:
         super().normalize(archive, logger)
 
@@ -223,12 +190,6 @@ class ImaginaryTime(Variables):
         """,
     )
 
-    def __init__(
-        self, m_def: 'Section' = None, m_context: 'Context' = None, **kwargs
-    ) -> None:
-        super().__init__(m_def, m_context, **kwargs)
-        self.name = self.m_def.name
-
     def normalize(self, archive: 'EntryArchive', logger: 'BoundLogger') -> None:
         super().normalize(archive, logger)
 
@@ -246,12 +207,6 @@ class KMesh(Variables):
         """,
     )
 
-    def __init__(
-        self, m_def: 'Section' = None, m_context: 'Context' = None, **kwargs
-    ) -> None:
-        super().__init__(m_def, m_context, **kwargs)
-        self.name = self.m_def.name
-
     def normalize(self, archive: 'EntryArchive', logger: 'BoundLogger') -> None:
         super().normalize(archive, logger)
 
@@ -266,11 +221,5 @@ class KLinePath(Variables):
         """,
     )
 
-    def __init__(
-        self, m_def: 'Section' = None, m_context: 'Context' = None, **kwargs
-    ) -> None:
-        super().__init__(m_def, m_context, **kwargs)
-        self.name = self.m_def.name
-
     def normalize(self, archive: 'EntryArchive', logger: 'BoundLogger') -> None:
         super().normalize(archive, logger)
diff --git a/tests/properties/test_band_gap.py b/tests/properties/test_band_gap.py
index df16d6bf4..bf8a63dc8 100644
--- a/tests/properties/test_band_gap.py
+++ b/tests/properties/test_band_gap.py
@@ -26,7 +26,6 @@ def test_default_quantities(self):
             electronic_band_gap.iri
             == 'http://fairmat-nfdi.eu/taxonomy/ElectronicBandGap'
         )
-        assert electronic_band_gap.name == 'ElectronicBandGap'
 
     @pytest.mark.parametrize(
         'momentum_transfer, type, result',
diff --git a/tests/properties/test_band_structure.py b/tests/properties/test_band_structure.py
index 34d735451..0e4c42678 100644
--- a/tests/properties/test_band_structure.py
+++ b/tests/properties/test_band_structure.py
@@ -33,7 +33,6 @@ def test_default_quantities(self, n_bands: Optional[int]):
             electronic_eigenvalues.iri
             == 'http://fairmat-nfdi.eu/taxonomy/ElectronicEigenvalues'
         )
-        assert electronic_eigenvalues.name == 'ElectronicEigenvalues'
 
     # @pytest.mark.parametrize(
     #     'occupation, result',
@@ -63,19 +62,19 @@ def test_default_quantities(self, n_bands: Optional[int]):
             (
                 None,
                 [[3, -2], [3, 1], [4, -2], [5, -1], [4, 0], [2, 0], [2, 1], [4, -3]],
-                False,
+                (),
                 (None, None),
             ),
             (
                 [[2, 2], [0, 0]],
                 [[3, -2], [3, 1], [4, -2], [5, -1], [4, 0], [2, 0], [2, 1], [4, -3]],
-                False,
+                (),
                 (None, None),
             ),  # `value` and `occupation` must have same shape
             (
                 [[0, 2], [0, 1], [0, 2], [0, 2], [0, 1.5], [0, 1.5], [0, 1], [0, 2]],
                 None,
-                False,
+                (),
                 (None, None),
             ),
             (
@@ -139,7 +138,7 @@ def test_order_eigenvalues(
         )
         order_result = electronic_eigenvalues.order_eigenvalues()
         if not order_result:
-            assert order_result == result_validation
+            assert result_validation == ()  # Empty tuple means validation failed
         else:
             sorted_value, sorted_occupation = order_result
             assert electronic_eigenvalues.m_cache['sorted_eigenvalues']
diff --git a/tests/properties/test_energies.py b/tests/properties/test_energies.py
deleted file mode 100644
index 83c7aa3f9..000000000
--- a/tests/properties/test_energies.py
+++ /dev/null
@@ -1,66 +0,0 @@
-from nomad_simulations.schema_packages.properties import (
-    FermiLevel,
-    KineticEnergy,
-    PotentialEnergy,
-    TotalEnergy,
-)
-
-
-class TestFermiLevel:
-    """
-    Test the `FermiLevel` class defined in `properties/energies.py`.
-    """
-
-    # ! Include this initial `test_default_quantities` method when testing your PhysicalProperty classes
-    def test_default_quantities(self):
-        """
-        Test the default quantities assigned when creating an instance of the `FermiLevel` class.
-        """
-        fermi_level = FermiLevel()
-        assert fermi_level.iri == 'http://fairmat-nfdi.eu/taxonomy/FermiLevel'
-        assert fermi_level.name == 'FermiLevel'
-
-
-class TestTotalEnergy:
-    """
-    Test the `TotalEnergy` class defined in `properties/energies.py`.
-    """
-
-    # ! Include this initial `test_default_quantities` method when testing your PhysicalProperty classes
-    def test_default_quantities(self):
-        """
-        Test the default quantities assigned when creating an instance of the `TotalEnergy` class.
-        """
-        total_energy = TotalEnergy()
-        # assert total_energy.iri == 'http://fairmat-nfdi.eu/taxonomy/TotalEnergy'
-        assert total_energy.name == 'TotalEnergy'
-
-
-class TestKineticEnergy:
-    """
-    Test the `KineticEnergy` class defined in `properties/energies.py`.
-    """
-
-    # ! Include this initial `test_default_quantities` method when testing your PhysicalProperty classes
-    def test_default_quantities(self):
-        """
-        Test the default quantities assigned when creating an instance of the `KineticEnergy` class.
-        """
-        kinetic_energy = KineticEnergy()
-        # assert kinetic_energy.iri == 'http://fairmat-nfdi.eu/taxonomy/KineticEnergy'
-        assert kinetic_energy.name == 'KineticEnergy'
-
-
-class TestPotentialEnergy:
-    """
-    Test the `PotentialEnergy` class defined in `properties/energies.py`.
-    """
-
-    # ! Include this initial `test_default_quantities` method when testing your PhysicalProperty classes
-    def test_default_quantities(self):
-        """
-        Test the default quantities assigned when creating an instance of the `PotentialEnergy` class.
-        """
-        potential_energy = PotentialEnergy()
-        # assert potential_energy.iri == 'http://fairmat-nfdi.eu/taxonomy/PotentialEnergy'
-        assert potential_energy.name == 'PotentialEnergy'
diff --git a/tests/properties/test_fermi_surface.py b/tests/properties/test_fermi_surface.py
index 78c105ebe..ebb5e0f88 100644
--- a/tests/properties/test_fermi_surface.py
+++ b/tests/properties/test_fermi_surface.py
@@ -24,4 +24,3 @@ def test_default_quantities(self, n_bands: Optional[int]):
         """
         fermi_surface = FermiSurface(n_bands=n_bands)
         assert fermi_surface.iri == 'http://fairmat-nfdi.eu/taxonomy/FermiSurface'
-        assert fermi_surface.name == 'FermiSurface'
diff --git a/tests/properties/test_forces.py b/tests/properties/test_forces.py
deleted file mode 100644
index 222fdf8e3..000000000
--- a/tests/properties/test_forces.py
+++ /dev/null
@@ -1,16 +0,0 @@
-from nomad_simulations.schema_packages.properties import TotalForce
-
-
-class TestTotalForce:
-    """
-    Test the `TotalForce` class defined in `properties/forces.py`.
-    """
-
-    # ! Include this initial `test_default_quantities` method when testing your PhysicalProperty classes
-    def test_default_quantities(self):
-        """
-        Test the default quantities assigned when creating an instance of the `TotalForce` class.
-        """
-        total_force = TotalForce()
-        # assert total_force.iri == 'http://fairmat-nfdi.eu/taxonomy/TotalForce'
-        assert total_force.name == 'TotalForce'
diff --git a/tests/properties/test_hopping_matrix.py b/tests/properties/test_hopping_matrix.py
index 72a4e5132..41bbf0b85 100644
--- a/tests/properties/test_hopping_matrix.py
+++ b/tests/properties/test_hopping_matrix.py
@@ -27,7 +27,6 @@ def test_default_quantities(self, n_orbitals: Optional[int], rank: Optional[list
         """
         hopping_matrix = HoppingMatrix(n_orbitals=n_orbitals)
         assert hopping_matrix.iri == 'http://fairmat-nfdi.eu/taxonomy/HoppingMatrix'
-        assert hopping_matrix.name == 'HoppingMatrix'
 
 
 class TestCrystalFieldSplitting:
@@ -51,4 +50,3 @@ def test_default_quantities(self, n_orbitals: Optional[int], rank: Optional[list
         assert (
             crystal_field.iri == 'http://fairmat-nfdi.eu/taxonomy/CrystalFieldSplitting'
         )
-        assert crystal_field.name == 'CrystalFieldSplitting'
diff --git a/tests/properties/test_permittivity.py b/tests/properties/test_permittivity.py
index 9d6574dd0..134a7a9a6 100644
--- a/tests/properties/test_permittivity.py
+++ b/tests/properties/test_permittivity.py
@@ -22,7 +22,6 @@ def test_default_quantities(self):
         """
         permittivity = Permittivity()
         assert permittivity.iri == 'http://fairmat-nfdi.eu/taxonomy/Permittivity'
-        assert permittivity.name == 'Permittivity'
 
     @pytest.mark.parametrize(
         'kmesh_grid, frequency, result',
diff --git a/tests/properties/test_spectral_profile.py b/tests/properties/test_spectral_profile.py
index ba6f2c963..2f592e225 100644
--- a/tests/properties/test_spectral_profile.py
+++ b/tests/properties/test_spectral_profile.py
@@ -34,7 +34,6 @@ def test_default_quantities(self):
             electronic_dos.iri
             == 'http://fairmat-nfdi.eu/taxonomy/ElectronicDensityOfStates'
         )
-        assert electronic_dos.name == 'ElectronicDensityOfStates'
 
     def test_resolve_energies_origin(self):
         """
@@ -209,8 +208,7 @@ def test_default_quantities(self):
         Test the default quantities assigned when creating an instance of the `AbsorptionSpectrum` class.
         """
         absorption_spectrum = AbsorptionSpectrum()
-        assert absorption_spectrum.iri is None  # Add iri when available
-        assert absorption_spectrum.name == 'AbsorptionSpectrum'
+        assert absorption_spectrum.iri == ''  # IRI is empty string by default
 
 
 class TestXASSpectrum:
@@ -218,15 +216,6 @@ class TestXASSpectrum:
     Test the `XASSpectrum` class defined in `properties/spectral_profile.py`.
     """
 
-    # ! Include this initial `test_default_quantities` method when testing your PhysicalProperty classes
-    def test_default_quantities(self):
-        """
-        Test the default quantities assigned when creating an instance of the `XASSpectrum` class.
-        """
-        xas_spectrum = XASSpectrum()
-        assert xas_spectrum.iri is None  # Add iri when available
-        assert xas_spectrum.name == 'XASSpectrum'
-
     @pytest.mark.parametrize(
         'xanes_energies, exafs_energies, xas_values',
         [
diff --git a/tests/test_outputs.py b/tests/test_outputs.py
index d5b9ba098..5ca611f30 100644
--- a/tests/test_outputs.py
+++ b/tests/test_outputs.py
@@ -21,39 +21,6 @@ class TestOutputs:
     Test the `Outputs` class defined in `outputs.py`.
     """
 
-    def test_number_of_properties(self):  # TODO: remove this test
-        """
-        Test how many properties are defined under `Outputs` and its order. This test is done in order to control better
-        which properties are already defined and in which order to control their normalizations
-        """
-        outputs = Outputs()
-        assert len(outputs.m_def.all_sub_sections) == 22
-        defined_properties = [
-            'fermi_levels',
-            'chemical_potentials',
-            'crystal_field_splittings',
-            'hopping_matrices',
-            'electronic_eigenvalues',
-            'electronic_band_gaps',
-            'electronic_dos',
-            'fermi_surfaces',
-            'electronic_band_structures',
-            'occupancies',
-            'electronic_greens_functions',
-            'electronic_self_energies',
-            'hybridization_functions',
-            'quasiparticle_weights',
-            'permittivities',
-            'absorption_spectra',
-            'xas_spectra',
-            'total_energies',
-            'kinetic_energies',
-            'potential_energies',
-            'total_forces',
-            'temperatures',
-        ]
-        assert list(outputs.m_def.all_sub_sections.keys()) == defined_properties
-
     @pytest.mark.parametrize(
         'band_gaps, values, result_length, result',
         [
@@ -95,7 +62,7 @@ def test_extract_spin_polarized_properties(
         outputs = Outputs()
 
         for i, band_gap in enumerate(band_gaps):
-            band_gap.value = [values[i]]
+            band_gap.value = values[i]
             outputs.electronic_band_gaps.append(band_gap)
         gaps = outputs.extract_spin_polarized_property(
             property_name='electronic_band_gaps'
@@ -103,7 +70,7 @@ def test_extract_spin_polarized_properties(
         assert len(gaps) == result_length
         if len(result) > 0:
             for i, result_gap in enumerate(result):
-                result_gap.value = [values[i]]
+                result_gap.value = values[i]
                 # ? comparing the sections does not work
                 assert gaps[i].value == result_gap.value
         else:
@@ -297,7 +264,7 @@ def test_get_last_scf_steps_value(
         for i, scf_step in enumerate(scf_last_steps):
             property_section = getattr(scf_step, 'electronic_band_gaps')
             if property_section is not None and values is not None:
-                property_section[i_property].value = [values[i]]
+                property_section[i_property].value = values[i]
         scf_values = scf_outputs.get_last_scf_steps_value(
             scf_last_steps=scf_last_steps,
             property_name='electronic_band_gaps',
@@ -314,8 +281,6 @@ def test_get_last_scf_steps_value(
             (0, None, '', 0, False),
             # no `self_consistency_ref` section
             (5, None, '', 0, False),
-            # no property matching `property_name`
-            (5, None, 'fermi_levels', 0, False),
             # `i_property` is out of range
             (5, None, 'electronic_band_gaps', 2, False),
             # property is not converged
diff --git a/tests/test_physical_properties.py b/tests/test_physical_properties.py
index 822504fea..41db2f749 100644
--- a/tests/test_physical_properties.py
+++ b/tests/test_physical_properties.py
@@ -1,17 +1,14 @@
-from typing import Optional, Union
+from typing import Callable
 
 import numpy as np
 import pytest
 from nomad.datamodel import EntryArchive
+from nomad.datamodel.metainfo.plot import PlotlyFigure
 from nomad.metainfo import Quantity
-from nomad.units import ureg
+from plotly.graph_objects import Figure
 
-from nomad_simulations.schema_packages.physical_property import (
-    PhysicalProperty,
-    validate_quantity_wrt_value,
-)
+from nomad_simulations.schema_packages.physical_property import PhysicalProperty
 
-# from nomad_simulations.schema_packages.variables import Variables
 from . import logger
 
 
@@ -21,90 +18,191 @@ class DummyPhysicalProperty(PhysicalProperty):
         unit='eV',
         shape=['*', '*', '*', '*'],
         description="""
-        This value is defined in order to test the `__setattr__` method in `PhysicalProperty`.
+        This value is defined in order to test the functionality in `PhysicalProperty`.
         """,
     )
 
+    def plot(self, **kwargs) -> list[PlotlyFigure]:
+        """Test implementation of plot method."""
+        fig = Figure()
+        fig.add_scatter(x=[1, 2, 3], y=[1, 4, 2], name='test')
+        plotly_figure = PlotlyFigure(label='test', figure=fig.to_plotly_json())
+        return [plotly_figure]
+
 
 class TestPhysicalProperty:
     """
     Test the `PhysicalProperty` class defined in `physical_property.py`.
     """
 
-    def test_setattr_value(self):
-        """
-        Test the `__setattr__` method when setting the `value` quantity of a physical property.
-        """
-        physical_property = DummyPhysicalProperty(
-            source='simulation',
-            # variables=[Variables(n_points=4), Variables(n_points=10)],
-        )
-        # `physical_property.value` must have full_shape=[4, 10, 3, 3]
-        value = np.ones((4, 10, 3, 3)) * ureg.eV
-        # assert physical_property.full_shape == list(value.shape)
-        physical_property.value = value
-        assert np.all(physical_property.value == value)
-
     def test_is_derived(self):
         """
         Test the `normalize` and `_is_derived` methods.
         """
         # Testing a directly parsed physical property
-        not_derived_physical_property = PhysicalProperty(source='simulation')
+        not_derived_physical_property = PhysicalProperty()
         assert not_derived_physical_property._is_derived() is False
         not_derived_physical_property.normalize(EntryArchive(), logger)
         assert not_derived_physical_property.is_derived is False
         # Testing a derived physical property
         derived_physical_property = PhysicalProperty(
-            source='analysis',
             physical_property_ref=not_derived_physical_property,
         )
         assert derived_physical_property._is_derived() is True
         derived_physical_property.normalize(EntryArchive(), logger)
         assert derived_physical_property.is_derived is True
 
+    def test_normalization_flag(self):
+        """
+        Test that the normalization flag prevents duplicate normalization.
+        """
+        property_obj = DummyPhysicalProperty()
+
+        # First normalization
+        property_obj.normalize(EntryArchive(), logger)
+        assert property_obj.m_cache.get('_is_normalized', False) is True
 
-# testing `validate_quantity_wrt_value` decorator
-class ValidatingClass:
-    def __init__(self, value=None, occupation=None):
-        self.value = value
-        self.occupation = occupation
-
-    @validate_quantity_wrt_value('occupation')
-    def validate_occupation(self) -> Union[bool, np.ndarray]:
-        return self.occupation
-
-
-@pytest.mark.parametrize(
-    'value, occupation, result',
-    [
-        (None, None, False),  # Both value and occupation are None
-        (np.array([[1, 2], [3, 4]]), None, False),  # occupation is None
-        (None, np.array([[0.5, 1], [0, 0.5]]), False),  # value is None
-        (np.array([[1, 2], [3, 4]]), np.array([]), False),  # occupation is empty
-        (
-            np.array([[1, 2], [3, 4]]),
-            np.array([[0.5, 1]]),
-            False,
-        ),  # Shapes do not match
-        (
-            np.array([[1, 2], [3, 4]]),
-            np.array([[0.5, 1], [0, 0.5]]),
-            np.array([[0.5, 1], [0, 0.5]]),
-        ),  # Valid case (return `occupation`)
-    ],
-)
-def test_validate_quantity_wrt_value(
-    value: Optional[np.ndarray],
-    occupation: Optional[np.ndarray],
-    result: Union[bool, np.ndarray],
-):
-    """
-    Test the `validate_quantity_wrt_value` decorator.
-    """
-    obj = ValidatingClass(value=value, occupation=occupation)
-    validation = obj.validate_occupation()
-    if isinstance(validation, bool):
-        assert validation == result
-    else:
-        assert np.allclose(validation, result)
+        # Store original figures count
+        original_figures_count = len(property_obj.figures)
+
+        # Second normalization should not duplicate work
+        property_obj.normalize(EntryArchive(), logger)
+
+        # Should still be marked as normalized
+        assert property_obj.m_cache.get('_is_normalized', False) is True
+        # Should not have duplicated figures
+        assert len(property_obj.figures) == original_figures_count
+
+    def test_plotting_and_contributions(self):
+        """
+        Test plotting integration and contributions normalization.
+        """
+        # Test main property plotting
+        property_obj = DummyPhysicalProperty()
+        property_obj.normalize(EntryArchive(), logger)
+
+        assert len(property_obj.figures) > 0
+        assert isinstance(property_obj.figures[0], PlotlyFigure)
+
+        # Test contributions
+        main_property = DummyPhysicalProperty()
+        contribution = DummyPhysicalProperty(name='contribution')
+        main_property.contributions = [contribution]
+        main_property.normalize(EntryArchive(), logger)
+
+        # Both should be normalized
+        assert main_property.m_cache.get('_is_normalized', False) is True
+        assert contribution.m_cache.get('_is_normalized', False) is True
+
+    @pytest.mark.parametrize(
+        'instantiator, reference',
+        [
+            (PhysicalProperty, 'PhysicalProperty'),
+            (DummyPhysicalProperty, 'DummyPhysicalProperty'),
+        ],
+    )
+    def test_name_setting_during_normalization(
+        self, instantiator: Callable, reference: str
+    ):
+        """
+        Test that the name is set during normalization for PhysicalProperty.
+        """
+        property_obj = instantiator()
+        property_obj.normalize(EntryArchive(), logger)
+        assert property_obj.name == reference
+
+    @pytest.mark.parametrize(
+        'has_nested_contributions, log_ref',
+        [(True, True), (False, False)],
+    )
+    def test_contributions_validation(
+        self, caplog, has_nested_contributions: bool, log_ref: bool
+    ):
+        """
+        Test contributions validation during normalization.
+
+        Args:
+            has_nested_contributions: Whether to create nested contribution structure
+            log_ref: Whether validation error should be logged
+        """
+        main_property = DummyPhysicalProperty(name='main')
+
+        if has_nested_contributions:
+            nested_contribution = DummyPhysicalProperty(name='nested')
+            contribution_with_nested = DummyPhysicalProperty(name='parent_contribution')
+            contribution_with_nested.contributions = [nested_contribution]
+            main_property.contributions = [contribution_with_nested]
+        else:
+            contribution1 = DummyPhysicalProperty(name='contrib1')
+            contribution2 = DummyPhysicalProperty(name='contrib2')
+            main_property.contributions = [contribution1, contribution2]
+
+        with caplog.at_level('ERROR'):
+            main_property.normalize(EntryArchive(), logger)
+
+        has_nested_error = any(
+            'nested contributions' in record.message.lower()
+            for record in caplog.records
+        )
+        assert has_nested_error == log_ref
+
+        if log_ref:
+            assert any('Contribution 0' in record.message for record in caplog.records)
+
+    @pytest.mark.parametrize(
+        'set_contribution_type_on_main, set_contribution_type_on_contrib, log_ref',
+        [
+            (False, False, False),
+            (False, True, False),
+            (True, False, True),
+            (True, True, True),
+        ],
+    )
+    def test_contribution_type_validation(
+        self,
+        caplog,
+        set_contribution_type_on_main: bool,
+        set_contribution_type_on_contrib: bool,
+        log_ref: bool,
+    ):
+        """
+        Test contribution_type validation during normalization.
+
+        Args:
+            set_contribution_type_on_main: Whether to set contribution_type on main property
+            set_contribution_type_on_contrib: Whether to set contribution_type on contribution
+            log_ref: Whether validation error should be logged
+        """
+        main_property = DummyPhysicalProperty(name='main')
+        if set_contribution_type_on_main:
+            main_property.contribution_type = 'invalid_main_type'
+
+        contribution = DummyPhysicalProperty(name='contrib')
+        if set_contribution_type_on_contrib:
+            contribution.contribution_type = 'valid_contrib_type'
+
+        main_property.contributions = [contribution]
+
+        with caplog.at_level('ERROR'):
+            main_property.normalize(EntryArchive(), logger)
+
+        has_contrib_type_error = any(
+            'contribution_type set but is not a contribution' in record.message.lower()
+            for record in caplog.records
+        )
+        assert has_contrib_type_error == log_ref
+
+    def test_is_contribution_method(self):
+        """
+        Test the _is_contribution helper method.
+        """
+        main_property = DummyPhysicalProperty(name='main')
+        contribution = DummyPhysicalProperty(name='contrib')
+
+        assert not main_property._is_contribution()
+        assert not contribution._is_contribution()
+
+        main_property.contributions = [contribution]
+        main_property.normalize(EntryArchive(), logger)
+
+        assert not main_property._is_contribution()
diff --git a/tests/test_variables.py b/tests/test_variables.py
index f9fc33126..9547a6740 100644
--- a/tests/test_variables.py
+++ b/tests/test_variables.py
@@ -33,3 +33,11 @@ def test_normalize(self, n_points: int, points: list, result: int):
         assert variable.get_n_points(logger) == result
         variable.normalize(EntryArchive(), logger)
         assert variable.n_points == result
+
+    def test_name_setting_during_normalization(self):
+        """
+        Test that the name is set during normalization for Variables.
+        """
+        variable = Variables()
+        variable.normalize(EntryArchive(), logger)
+        assert variable.name == 'Variables'