PyMPDATA_examples.utils.discretisation 

 1import numpy as np
 2from scipy import integrate
 3
 4
 5def from_pdf_2d(pdf: callable, xrange: list, yrange: list, gridsize: list):
 6    z = np.empty(gridsize)
 7    dx, dy = (xrange[1] - xrange[0]) / gridsize[0], (yrange[1] - yrange[0]) / gridsize[
 8        1
 9    ]
10    for i in range(gridsize[0]):
11        for j in range(gridsize[1]):
12            z[i, j] = (
13                integrate.nquad(
14                    pdf,
15                    ranges=(
16                        (xrange[0] + dx * i, xrange[0] + dx * (i + 1)),
17                        (yrange[0] + dy * j, yrange[0] + dy * (j + 1)),
18                    ),
19                )[0]
20                / dx
21                / dy
22            )
23    x = np.linspace(xrange[0] + dx / 2, xrange[1] - dx / 2, gridsize[0])
24    y = np.linspace(yrange[0] + dy / 2, yrange[1] - dy / 2, gridsize[1])
25    return x, y, z
26
27
28def from_cdf_1d(cdf: callable, x_min: float, x_max: float, nx: int):
29    dx = (x_max - x_min) / nx
30    x = np.linspace(x_min + dx / 2, x_max - dx / 2, nx)
31    xh = np.linspace(x_min, x_max, nx + 1)
32    y = np.diff(cdf(xh)) / dx
33    return x, y
34
35
36def discretised_analytical_solution(rh, pdf_t, midpoint_value=False, r=None):
37    if midpoint_value:
38        assert r is not None
39    else:
40        assert r is None
41    output = np.empty(rh.shape[0] - 1)
42    for i in range(output.shape[0]):
43        if midpoint_value:
44            output[i] = pdf_t(r[i])
45        else:
46            dcdf, _ = integrate.quad(pdf_t, rh[i], rh[i + 1])
47            output[i] = dcdf / (rh[i + 1] - rh[i])
48    return output
def from_pdf_2d( pdf: <built-in function callable>, xrange: list, yrange: list, gridsize: list): 
 6def from_pdf_2d(pdf: callable, xrange: list, yrange: list, gridsize: list):
 7    z = np.empty(gridsize)
 8    dx, dy = (xrange[1] - xrange[0]) / gridsize[0], (yrange[1] - yrange[0]) / gridsize[
 9        1
10    ]
11    for i in range(gridsize[0]):
12        for j in range(gridsize[1]):
13            z[i, j] = (
14                integrate.nquad(
15                    pdf,
16                    ranges=(
17                        (xrange[0] + dx * i, xrange[0] + dx * (i + 1)),
18                        (yrange[0] + dy * j, yrange[0] + dy * (j + 1)),
19                    ),
20                )[0]
21                / dx
22                / dy
23            )
24    x = np.linspace(xrange[0] + dx / 2, xrange[1] - dx / 2, gridsize[0])
25    y = np.linspace(yrange[0] + dy / 2, yrange[1] - dy / 2, gridsize[1])
26    return x, y, z
def from_cdf_1d( cdf: <built-in function callable>, x_min: float, x_max: float, nx: int): 
29def from_cdf_1d(cdf: callable, x_min: float, x_max: float, nx: int):
30    dx = (x_max - x_min) / nx
31    x = np.linspace(x_min + dx / 2, x_max - dx / 2, nx)
32    xh = np.linspace(x_min, x_max, nx + 1)
33    y = np.diff(cdf(xh)) / dx
34    return x, y
def discretised_analytical_solution(rh, pdf_t, midpoint_value=False, r=None): 
37def discretised_analytical_solution(rh, pdf_t, midpoint_value=False, r=None):
38    if midpoint_value:
39        assert r is not None
40    else:
41        assert r is None
42    output = np.empty(rh.shape[0] - 1)
43    for i in range(output.shape[0]):
44        if midpoint_value:
45            output[i] = pdf_t(r[i])
46        else:
47            dcdf, _ = integrate.quad(pdf_t, rh[i], rh[i + 1])
48            output[i] = dcdf / (rh[i + 1] - rh[i])
49    return output