A method for predicting a state-of-power, SoP, of an electric energy storage system, ESS, comprising at least two battery units (2, 2', 2") electrically connected in parallel, the method comprising: - obtaining (S1-1) operational data from the at least two battery units of the ESS during operation of the ESS; - computing (S1-2) the state-of-power of the ESS based on the obtained operational data and by using an algorithm based on a system-level model of the ESS, wherein the system-level model of the ESS takes into account on one hand each one of the at least two battery units (2, 2', 2") of the ESS, and on the other hand at least one electrical connection between the at least two battery units (2, 2', 2"), and wherein the system-level model of the ESS further takes into account a dynamic parallel load distribution between the at least two battery units (2, 2', 2"), the method further characterized in that operating limits of at least one constrained variable of each one of the at least two battery units (2, 2', 2") are used as input to the system-level model, wherein the at least one constrained variable includes at least one of battery unit current, battery unit terminal voltage, battery unit temperature, battery unit state-of-charge, and battery unit open circuit voltage, the step of computing (S1-2) the state-of-power of the ESS comprises solving a constrained optimization problem, in which a possible load current magnitude and/or a possible load power magnitude for the ESS is/are maximised subject to the operating limits of the at least one constrained variable, the state-of-power is predicted for a predefinable prediction time horizon ([ t 0, t 0 + Δt]), and the estimation comprises predicting an evolution of the at least one constrained variable during the prediction time horizon ([ t 0, t 0 + Δ t]), a maximum possible load current magnitude and/or a maximum possible load power magnitude for the ESS is/are set to be constant over the prediction time horizon ([ t 0, t 0 + Δ t]), and for each individual battery unit (2, 2', 2"), a battery unit load power or load current is allowed to vary over time during the prediction time horizon ([ t 0, t 0 + Δ t]). - obtaining (S1-1) operational data from the at least two battery units of the ESS during operation of the ESS; - computing (S1-2) the state-of-power of the ESS based on the obtained operational data and by using an algorithm based on a system-level model of the ESS, wherein the system-level model of the ESS takes into account on one hand each one of the at least two battery units (2, 2', 2") of the ESS, and on the other hand at least one electrical connection between the at least two battery units (2, 2', 2"), and wherein the system-level model of the ESS further takes into account a dynamic parallel load distribution between the at least two battery units (2, 2', 2"), operating limits of at least one constrained variable of each one of the at least two battery units (2, 2', 2") are used as input to the system-level model, wherein the at least one constrained variable includes at least one of battery unit current, battery unit terminal voltage, battery unit temperature, battery unit state-of-charge, and battery unit open circuit voltage, the step of computing (S1-2) the state-of-power of the ESS comprises solving a constrained optimization problem, in which a possible load current magnitude and/or a possible load power magnitude for the ESS is/are maximised subject to the operating limits of the at least one constrained variable, the state-of-power is predicted for a predefinable prediction time horizon ([ t 0, t 0 + Δt]), and the estimation comprises predicting an evolution of the at least one constrained variable during the prediction time horizon ([ t 0, t 0 + Δ t]), a maximum possible load current magnitude and/or a maximum possible load power magnitude for the ESS is/are set to be constant over the prediction time horizon ([ t 0, t 0 + Δ t]), and for each individual battery unit (2, 2', 2"), a battery unit load power or load current is allowed to vary over time during the prediction time horizon ([ t 0, t 0 + Δ t]).