Skip to content

Commit

Permalink
ruff
Browse files Browse the repository at this point in the history
  • Loading branch information
cneyens committed Jan 18, 2025
1 parent 92af682 commit f3e5fc8
Show file tree
Hide file tree
Showing 11 changed files with 557 additions and 305 deletions.
10 changes: 10 additions & 0 deletions .github/workflows/gha-ruff.yml
Original file line number Diff line number Diff line change
@@ -0,0 +1,10 @@
name: gha-ruff
on: [ push, pull_request ]
jobs:
ruff:
runs-on: ubuntu-latest
steps:
- uses: actions/checkout@v4
- uses: astral-sh/ruff-action@v3
- run: ruff check
- run: ruff format --check
3 changes: 1 addition & 2 deletions adepy/__init__.py
Original file line number Diff line number Diff line change
@@ -1,7 +1,6 @@
# ruff : noqa: F401
from ._version import __version__

from adepy.oneD import finite1, finite3, seminf1, seminf3
from adepy.twoD import point2, stripf, stripi, gauss
from adepy.threeD import point3, patchi, patchf


15 changes: 8 additions & 7 deletions adepy/_helpers.py
Original file line number Diff line number Diff line change
Expand Up @@ -9,25 +9,27 @@
functype = ctypes.CFUNCTYPE(ctypes.c_double, ctypes.c_double)
erfc_fn = functype(addr)

@vectorize('float64(float64)')

@vectorize("float64(float64)")
def _vec_erfc(x):
return erfc_fn(x)


@njit
def _erfc_nb(x):
return _vec_erfc(x)

def _integrate(integrand, t, *args, order=100, method='legendre'):

if method == 'legendre':
def _integrate(integrand, t, *args, order=100, method="legendre"):
if method == "legendre":
roots, weights = roots_legendre(order)

def integrate(t, *args):
roots_adj = roots * (t - 0) / 2 + (0 + t) / 2
F = integrand(roots_adj, *args).dot(weights) * (t - 0)/2
F = integrand(roots_adj, *args).dot(weights) * (t - 0) / 2
return F
elif method == 'quadrature':

elif method == "quadrature":

def integrate(t, *args):
F = quad(integrand, 0, t, args=args.items)
Expand All @@ -39,4 +41,3 @@ def integrate(t, *args):
integrate_vec = np.vectorize(integrate)
term = integrate_vec(t, *args)
return term

2 changes: 1 addition & 1 deletion adepy/_version.py
Original file line number Diff line number Diff line change
@@ -1 +1 @@
__version__ = "0.1.0"
__version__ = "0.1.0"
106 changes: 73 additions & 33 deletions adepy/oneD.py
Original file line number Diff line number Diff line change
@@ -1,41 +1,55 @@
from scipy.optimize import brentq
import numpy as np
from numba import njit, vectorize
from adepy._helpers import _erfc_nb as erfc


# @njit
def _bserie_finite1(betas, x, t, Pe, L, D, lamb):
bs = 0.0

if lamb == 0:
for b in betas:
bs += b * np.sin(b * x / L) * (b**2 + (Pe / 2)**2) * np.exp(-b**2 * D * t / L**2) /\
((b**2 + (Pe / 2)**2 + Pe / 2) * (b**2 + (Pe / 2)**2 + lamb / D * L**2))
bs += (
b
* np.sin(b * x / L)
* (b**2 + (Pe / 2) ** 2)
* np.exp(-(b**2) * D * t / L**2)
/ (
(b**2 + (Pe / 2) ** 2 + Pe / 2)
* (b**2 + (Pe / 2) ** 2 + lamb / D * L**2)
)
)
else:
for b in betas:
bs += b * np.sin(b * x / L) * np.exp(-b**2 * D * t / L**2) / (b**2 + (Pe / 2)**2 + Pe / 2)
bs += (
b
* np.sin(b * x / L)
* np.exp(-(b**2) * D * t / L**2)
/ (b**2 + (Pe / 2) ** 2 + Pe / 2)
)
return bs


def finite1(c0, x, t, v, al, L, Dm=0, lamb=0, R=1.0, nterm=1000):

x = np.atleast_1d(x)
t = np.atleast_1d(t)

if len(x) > 1 and len(t) > 1:
raise ValueError('Either x or t should have length 1')
raise ValueError("Either x or t should have length 1")

D = al * v + Dm

# apply retardation coefficient to right-hand side
v = v / R
D = D / R

Pe = v * L / D

# find roots
def betaf(b):
return b * 1 / np.tan(b) + Pe / 2

intervals = [np.pi * i for i in range(nterm)]
intervals = [np.pi * i for i in range(nterm)]
betas = []
for i in range(len(intervals) - 1):
mi = intervals[i] + 1e-10
Expand All @@ -50,38 +64,47 @@ def betaf(b):
term0 = 1.0
else:
u = np.sqrt(v**2 + 4 * lamb * D)
term0 = (np.exp((v - u) * x / (2 * D)) + (u - v) / (u + v) * np.exp((v + u) * x / (2 * D) - u * L / D)) /\
(1 + (u - v) / (u + v) * np.exp(-u * L / D))

term0 = (
np.exp((v - u) * x / (2 * D))
+ (u - v) / (u + v) * np.exp((v + u) * x / (2 * D) - u * L / D)
) / (1 + (u - v) / (u + v) * np.exp(-u * L / D))

term1 = -2 * np.exp(v * x / (2 * D) - v**2 * t / (4 * D) - lamb * t)

return c0 * (term0 + term1 * series)


# @njit
def _bserie_finite3(betas, x, t, Pe, L, D, lamb):
bs = 0.0
for b in betas:
bs += b * (b * np.cos(b * x / L) + (Pe / 2) * np.sin(b * x / L)) / (b**2 + (Pe / 2)**2 + Pe) *\
np.exp(-b**2 * D * t / L**2) / (b**2 + (Pe / 2)**2 + lamb * L**2 / D)
bs += (
b
* (b * np.cos(b * x / L) + (Pe / 2) * np.sin(b * x / L))
/ (b**2 + (Pe / 2) ** 2 + Pe)
* np.exp(-(b**2) * D * t / L**2)
/ (b**2 + (Pe / 2) ** 2 + lamb * L**2 / D)
)
return bs


def finite3(c0, x, t, v, al, L, Dm=0, lamb=0, R=1.0, nterm=1000):
# https://github.com/BYL4746/columntracer/blob/main/columntracer.py

x = np.atleast_1d(x)
t = np.atleast_1d(t)

if len(x) > 1 and len(t) > 1:
raise ValueError('Either x or t should have length 1')
raise ValueError("Either x or t should have length 1")

D = al * v + Dm

# apply retardation coefficient to right-hand side
v = v / R
D = D / R

Pe = v * L / D

# find roots
def betaf(b):
return b * 1 / np.tan(b) - b**2 / Pe + Pe / 4
Expand All @@ -101,13 +124,16 @@ def betaf(b):
term0 = 1.0
else:
u = np.sqrt(v**2 + 4 * lamb * D)
term0 = (np.exp((v - u) * x / (2 * D)) + (u - v) / (u + v) * np.exp((v + u) * x / (2 * D) - u * L / D)) /\
((u + v) / (2 * v) - (u - v)**2 / (2 * v * (u + v)) * np.exp(-u * L / D))

term0 = (
np.exp((v - u) * x / (2 * D))
+ (u - v) / (u + v) * np.exp((v + u) * x / (2 * D) - u * L / D)
) / ((u + v) / (2 * v) - (u - v) ** 2 / (2 * v * (u + v)) * np.exp(-u * L / D))

term1 = -2 * Pe * np.exp(v * x / (2 * D) - v**2 * t / (4 * D) - lamb * t)

return c0 * (term0 + term1 * series)


def seminf1(c0, x, t, v, al, Dm=0, lamb=0, R=1.0):
x = np.atleast_1d(x)
t = np.atleast_1d(t)
Expand All @@ -119,11 +145,13 @@ def seminf1(c0, x, t, v, al, Dm=0, lamb=0, R=1.0):
D = D / R

u = np.sqrt(v**2 + 4 * lamb * D)
term = np.exp(x * (v - u) / (2 * D)) * erfc((x - u * t) / (2 * np.sqrt(D * t))) + \
np.exp(x * (v + u) / (2 * D)) * erfc((x + u * t) / (2 * np.sqrt(D * t)))
term = np.exp(x * (v - u) / (2 * D)) * erfc(
(x - u * t) / (2 * np.sqrt(D * t))
) + np.exp(x * (v + u) / (2 * D)) * erfc((x + u * t) / (2 * np.sqrt(D * t)))

return c0 * 0.5 * term


def seminf3(c0, x, t, v, al, Dm=0, lamb=0, R=1.0):
x = np.atleast_1d(x)
t = np.atleast_1d(t)
Expand All @@ -133,17 +161,29 @@ def seminf3(c0, x, t, v, al, Dm=0, lamb=0, R=1.0):
# apply retardation coefficient to right-hand side
v = v / R
D = D / R

u = np.sqrt(v**2 + 4 * lamb * D)
if lamb == 0:
term = 0.5 * erfc((x - v * t) / (2 * np.sqrt(D * t))) + np.sqrt(t * v**2 / (np.pi * D)) * np.exp(-(x - v * t)**2 / (4 * D * t)) - \
0.5 * (1 + v * x / D + t * v**2 / D) * np.exp(v * x / D) * erfc((x + v * t) / (2 * np.sqrt(D * t)))
term = (
0.5 * erfc((x - v * t) / (2 * np.sqrt(D * t)))
+ np.sqrt(t * v**2 / (np.pi * D))
* np.exp(-((x - v * t) ** 2) / (4 * D * t))
- 0.5
* (1 + v * x / D + t * v**2 / D)
* np.exp(v * x / D)
* erfc((x + v * t) / (2 * np.sqrt(D * t)))
)
term0 = 1.0
else:
term = 2 * np.exp(x * v / D - lamb * t) * erfc((x + v * t) / (2 * np.sqrt(D * t))) +\
(u / v - 1) * np.exp(x * (v - u) / (2 * D)) * erfc((x - u * t) / (2 * np.sqrt(D * t))) -\
(u / v - 1) * np.exp(x * (v + u) / (2 * D)) * erfc((x + u * t) / (2 * np.sqrt(D * t)))
term = (
2 * np.exp(x * v / D - lamb * t) * erfc((x + v * t) / (2 * np.sqrt(D * t)))
+ (u / v - 1)
* np.exp(x * (v - u) / (2 * D))
* erfc((x - u * t) / (2 * np.sqrt(D * t)))
- (u / v - 1)
* np.exp(x * (v + u) / (2 * D))
* erfc((x + u * t) / (2 * np.sqrt(D * t)))
)
term0 = v**2 / (4 * lamb * D)

return c0 * term0 * term

Loading

0 comments on commit f3e5fc8

Please sign in to comment.