Vectorization

[juan11iguel]

Hello, first of all I would like to thank you for the great library, I use it in all my python projects.

The main limitation that I find in the library is that very often we work with arrays, and whenever a water property needs to be computed, having to loop trough each element to get it is a performance bottle neck. I wonder if there could be any possibility of implementing a way of giving as input a vector and for the library to return a vector of water properties (implemented with vectorization internally)

Best regards, Juanmi

[jjgomera]

Hi Juanmi, That's would be a great enhancement, I give a try next week

[juan11iguel]

That's good to hear!

Though after a bit of reading and digging in the library I am not sure it's a simple task. As I understand it, internally a solver is used, and I am not sure how that could be vectorized..

Also, maybe you find this source useful:

There are quite a few implementations of IAPWS-95 and IAPWS-97 out there. Besides the many commercial implementations, the are the following excellent open source ones:

• iapws by Juan José Gómez Romera, GPL3 licensed, containing IAPWS-95 and IAPWS-97 among other standards. Implemented in Python.
• CoolProp by Ian Bell, MIT licensed and containing IAPWS-95 and IAPWS-97 along with their transport properties. Implemented in C++ with an excellent interface to Python among other languages.
• freesteam by John Pye, GPL3 licensed, containing most of IAPWS-97 and the transport properties. Implemented in C.

There are many more, but these are the best developed libraries that can be used from Python. Water is so common and present in so many calculations that for many applications it is important to make it as fast as possible. IAPWS-95 is conventionally slow; properties are requested at a specified temperature T and pressure P, but the equation of state’s input variables are temperature and density! A numerical solver must be used in this case to find the density which yields the specified pressure. This density-solution procedure is normally the slowest part, although computing some properties requires many derivatives that can be slow also. A good conventional density solver will take ~10-30 μs on a modern computer. Only the CPU clockspeed really matters for this calculation time. It was discovered that with the use of Common subexpression elimination, the calculation could be speed up quite a lot. Additionally, if the IAPWS-95 density solution is initialized by the IAPWS-97 explicit calculation (applicable most of the time but not always), a few more iterations can be saved. The net result of these optimizations is a greatly improved density solve time - normally 2.5-4 μs when running with PyPy or Numba. The con to this approach is that the code is nearly unreadable, and it would not be possible to update the coefficients without rewriting the implementation. As IAPWS-95 is a static model which will be the best one available for many years to come, this is an acceptable trade off.

Source

[jjgomera]

Added multithreading support for IAPWS95 (valid for D2O too), very useful in isoline calculations, in my equipment i have a 6x speed up.

    from iapws import IAPWS95
    from numpy import arange
    from time import time

    x = arange(0, 1.01, 0.01)

    def fi(x):
        return IAPWS95(P=20.8, x=x)

    start = time()
    for xi in x:
        fi(xi)
    print(f'Without multiprocessing: {time() - start}')

    start = time()
    states = IAPWS95.from_list("P", 20.8, "x", x)
    print(f'With multiprocessing: {time() - start}')

For IAPWS97 it's necessary many more code rewriting i think.

[juan11iguel]

Hi, sorry for the late reply, I tried the latest commit and I can confirm everything works as expected and getting about the same performance increases.

I tested it the following way:

1. Install uv to manage dependencies (curl -LsSf https://astral.sh/uv/install.sh | sh)
2. Run the script test_iapws.py with uv run test_iapws.py

Reading inline script metadata from: test_iapws.py
 Updated https://github.com/jjgomera/iapws.git (c9d2767)
   Built iapws @ git+https://github.com/jjgomera/iapws.git@c9d2767f03a84d067ab17af6fd9ed2cc5322169d
Installed 3 packages in 111ms

Average improvement of using multiprocessing after 6 runs: 6.2X speed increase

test_iapws.py content:

# /// script
# dependencies = [
#   "numpy",
#   "iapws @ git+https://github.com/jjgomera/iapws.git@c9d2767f03a84d067ab17af6fd9ed2cc5322169d",
# ]
# ///

from iapws import IAPWS95
import numpy as np
from time import time

def fi(x):
    return IAPWS95(P=20.8, x=x)

def evaluate(x: np.ndarray[float], n_evals: int = 3, print_result: bool = True) -> list[tuple[float, float]]:

    elapsed_times: list[tuple[float, float]] = []

    for _ in range(n_evals):
        start = time()
        states = [fi(xi) for xi in x]
        elapsed_time = time() - start

        start = time()
        states = IAPWS95.from_list("P", 20.8, "x", x)
        elapsed_time_parallel = time() - start

        elapsed_times.append( (elapsed_time, elapsed_time_parallel) )

        if print_result:
            print(f'Without multiprocessing: {elapsed_time:.2f} seconds')
            print(f'With multiprocessing: {elapsed_time_parallel:.2f} seconds, improvement of {elapsed_time/elapsed_time_parallel:.1f}X')

        return elapsed_times

def main() -> None:
    x = np.arange(0, 1.01, 0.01)
    n_evals = 10
    print_result = False    

    elapsed_times = evaluate(x, n_evals=n_evals, print_result=print_result)
    avg_improvement = sum(et / etp for et, etp in elapsed_times) / len(elapsed_times)

    print(f"Average improvement of using multiprocessing after {n_evals} runs: {avg_improvement:.1f}X speed increase")

if __name__ == "__main__":
    main()

Gracias!