PyMPDATA_examples.Olesik_et_al_2022.equilibrium_drop_growth 

 1import numpy as np
 2
 3
 4class DrDt:
 5    """eq. 7.20 in Rogers and Yau 1989"""
 6
 7    def __init__(self, ksi_1, S):
 8        self.ksi = (S - 1) * ksi_1
 9
10    def __call__(self, r):
11        return self.ksi / r
12
13    def mean(self, r1, r2):
14        return self.ksi * np.log(r2 / r1) / (r2 - r1)
15
16
17# pylint: disable=too-few-public-methods
18class PdfEvolver:
19    """eq. 7.32 in Rogers and Yau 1989"""
20
21    def __init__(self, pdf, drdt: DrDt, t):
22        self.t = t
23        self.pdf = pdf
24        self.drdt = drdt
25
26    def __call__(self, r):
27        with np.errstate(invalid="ignore"):
28            arg = np.sqrt(r**2 - 2 * self.drdt.ksi * self.t)
29        result = r / arg * self.pdf(arg)
30
31        if isinstance(result.magnitude, np.ndarray):
32            result = (
33                np.where(np.isfinite(result.magnitude), result.magnitude, 0)
34                * result.units
35            )
36        else:
37            if not np.isfinite(result):
38                result = 0 * result.units
39
40        return result
class DrDt: 
 5class DrDt:
 6    """eq. 7.20 in Rogers and Yau 1989"""
 7
 8    def __init__(self, ksi_1, S):
 9        self.ksi = (S - 1) * ksi_1
10
11    def __call__(self, r):
12        return self.ksi / r
13
14    def mean(self, r1, r2):
15        return self.ksi * np.log(r2 / r1) / (r2 - r1)

eq. 7.20 in Rogers and Yau 1989

DrDt(ksi_1, S) 
8    def __init__(self, ksi_1, S):
9        self.ksi = (S - 1) * ksi_1
ksi
def mean(self, r1, r2): 
14    def mean(self, r1, r2):
15        return self.ksi * np.log(r2 / r1) / (r2 - r1)
class PdfEvolver: 
19class PdfEvolver:
20    """eq. 7.32 in Rogers and Yau 1989"""
21
22    def __init__(self, pdf, drdt: DrDt, t):
23        self.t = t
24        self.pdf = pdf
25        self.drdt = drdt
26
27    def __call__(self, r):
28        with np.errstate(invalid="ignore"):
29            arg = np.sqrt(r**2 - 2 * self.drdt.ksi * self.t)
30        result = r / arg * self.pdf(arg)
31
32        if isinstance(result.magnitude, np.ndarray):
33            result = (
34                np.where(np.isfinite(result.magnitude), result.magnitude, 0)
35                * result.units
36            )
37        else:
38            if not np.isfinite(result):
39                result = 0 * result.units
40
41        return result

eq. 7.32 in Rogers and Yau 1989

PdfEvolver( pdf, drdt: DrDt, t) 
22    def __init__(self, pdf, drdt: DrDt, t):
23        self.t = t
24        self.pdf = pdf
25        self.drdt = drdt
t
pdf
drdt