diff --git a/.gitlab-ci.yml b/.gitlab-ci.yml index 8de78af3089257c96dcb89f1d1d543f1234446c8..906e8e4b84476fa82cca82f8ed4c6ab010137c24 100644 --- a/.gitlab-ci.yml +++ b/.gitlab-ci.yml @@ -10,6 +10,7 @@ Python 3.5 POCL: - pocl except: - tags + Documentation: script: - EXTRA_INSTALL="numpy mako" @@ -19,3 +20,12 @@ Documentation: - python3.5 only: - master + +Flake8: + script: + - curl -L -O -k https://gitlab.tiker.net/inducer/ci-support/raw/master/prepare-and-run-flake8.sh + - ". ./prepare-and-run-flake8.sh pytential test" + tags: + - python3.5 + except: + - tags diff --git a/pytential/symbolic/old_diffop_primitives.py b/pytential/symbolic/old_diffop_primitives.py new file mode 100644 index 0000000000000000000000000000000000000000..23028e14f635a6eb44b038742637578c00cc6d1c --- /dev/null +++ b/pytential/symbolic/old_diffop_primitives.py @@ -0,0 +1,173 @@ + +# {{{ differential operators on layer potentials + +def grad_S(kernel, arg, dim): + from pytools.obj_array import log_shape + arg_shape = log_shape(arg) + result = np.zeros(arg_shape+(dim,), dtype=object) + from pytools import indices_in_shape + for i in indices_in_shape(arg_shape): + for j in range(dim): + result[i+(j,)] = IntGdTarget(kernel, arg[i], j) + return result + + +def grad_D(kernel, arg, dim): + from pytools.obj_array import log_shape + arg_shape = log_shape(arg) + result = np.zeros(arg_shape+(dim,), dtype=object) + from pytools import indices_in_shape + for i in indices_in_shape(arg_shape): + for j in range(dim): + result[i+(j,)] = IntGdMixed(kernel, arg[i], j) + return result + + +def tangential_surf_grad_source_S(kernel, arg, dim=3): + from pytools.obj_array import make_obj_array + return make_obj_array([ + IntGdSource(kernel, arg, + ds_direction=make_tangent(i, dim, "src")) + for i in range(dim-1)]) + + +def surf_grad_S(kernel, arg, dim): + """ + :arg dim: The dimension of the ambient space. + """ + + return project_to_tangential(cse(grad_S(kernel, arg, dim))) + + +def div_S_volume(kernel, arg): + return sum(IntGdTarget(kernel, arg_n, n) for n, arg_n in enumerate(arg)) + + +def curl_S_volume(kernel, arg): + from pytools import levi_civita + from pytools.obj_array import make_obj_array + + return make_obj_array([ + sum( + levi_civita((l, m, n)) * IntGdTarget(kernel, arg[n], m) + for m in range(3) for n in range(3)) + for l in range(3)]) + + +def curl_curl_S_volume(k, arg): + # By vector identity, this is grad div S volume + k^2 S_k(arg), + # since S_k(arg) satisfies a Helmholtz equation. + + from pytools.obj_array import make_obj_array + + def swap_min_first(i, j): + if i < j: + return i, j + else: + return j, i + + return make_obj_array([ + sum(IntGd2Target(k, arg[m], *swap_min_first(m, n)) for m in range(3)) + for n in range(3)]) + k**2*S(k, arg) + + +def nxcurl_S(kernel, loc, arg): + """ + :arg loc: one of three values: + * +1 on the side of the surface toward which + the normal points ('exterior' of the surface), + * 0 on the surface, or evaluated at a volume target + * -1 on the interior of the surface. + """ + nxcurl_S = np.cross(normal(3), curl_S_volume(kernel, arg)) + assert loc in [-1, 0, 1], "invalid value for 'loc' (%s)" % loc + return nxcurl_S + loc*(1/2)*arg + + +def surface_laplacian_S_squared(u, invertibility_scale=0): + """ + :arg u: The field to which the surface Laplacian is applied. + """ + # http://wiki.tiker.net/HellsKitchen/SurfaceLaplacian + + Su = cse(S(0, u), "su_from_surflap") + + return ( + - 2*mean_curvature()*Sp(0, Su) + - ((Spp(0, Su)+Dp(0, Su))-(-1/4*u+Sp(0, Sp(0, u)))) + - invertibility_scale * mean(S(0, Su))*Ones()) + + +def S_surface_laplacian_S(u, dim, invertibility_scale=0, qbx_fix_scale=0): + """ + :arg u: The field to which the surface Laplacian is applied. + """ + + # This converges, but appears to work quite poorly compared to the above. + + tgrad_Su = cse( + project_to_tangential(grad_S(0, u, dim)), + "tgrad_su_from_surflap") + + return ( + - IntGdSource(0, Ones(), ds_direction=real(tgrad_Su)) + - 1j*IntGdSource(0, Ones(), ds_direction=imag(tgrad_Su)) + - invertibility_scale * S(0, Ones()*mean(S(0, u))) + - qbx_fix_scale * ( + u + # D+ - D- = identity (but QBX will compute the + # 'compact part of the identity' -- call that I*) + - ( + D(0, u, qbx_forced_limit=+1) + - D(0, u, qbx_forced_limit=-1)) + ) + # The above is I - I*, which means only the high-frequency + # bits of the identity are left. + ) + +# }}} + + +# {{{ geometric operations + +def xyz_to_tangential(xyz_vec, which=None): + d = len(xyz_vec) + x2l = xyz_to_local_matrix(d) + return np.dot(x2l[:-1], xyz_vec) + + +def tangential_to_xyz(tangential_vec, which=None): + d = len(tangential_vec) + 1 + x2l = xyz_to_local_matrix(d) + return np.dot(x2l[:-1].T, tangential_vec) + + +def project_to_tangential(xyz_vec, which=None): + return tangential_to_xyz( + cse(xyz_to_tangential(xyz_vec, which), which)) + + +def n_dot(vec, which=None): + return np.dot(normal(len(vec), which), vec) + + +def n_cross(vec, which=None): + nrm = normal(3, which) + + from pytools import levi_civita + from pytools.obj_array import make_obj_array + return make_obj_array([ + sum( + levi_civita((i, j, k)) * nrm[j] * vec[k] + for j in range(3) for k in range(3)) + for i in range(3)]) + + +def surf_n_cross(tangential_vec): + assert len(tangential_vec) == 2 + from pytools.obj_array import make_obj_array + return make_obj_array([-tangential_vec[1], tangential_vec[0]]) + +# }}} + + diff --git a/pytential/symbolic/pde/scalar.py b/pytential/symbolic/pde/scalar.py index 94cc39a9923d8c7cacb87d9ca03828a6b613480a..8cb2f825bc11b9449e02a682bddc5661edfadb3c 100644 --- a/pytential/symbolic/pde/scalar.py +++ b/pytential/symbolic/pde/scalar.py @@ -227,13 +227,13 @@ class NeumannOperator(L2WeightedPDEOperator): cse(S(lknl, inv_sqrt_w_u))) if self.use_improved_operator: - Dp0S0u = -0.25*u + Sp( + Dp0S0u = -0.25*u + Sp( # noqa lknl, # noqa Sp(lknl, inv_sqrt_w_u, qbx_forced_limit="avg"), qbx_forced_limit="avg") if isinstance(self.kernel, HelmholtzKernel): - DpS0u = ( + DpS0u = ( # noqa Dp(self.kernel_and_args - lknl, # noqa cse(S(lknl, inv_sqrt_w_u, qbx_forced_limit=+1)), qbx_forced_limit=+1) diff --git a/pytential/symbolic/primitives.py b/pytential/symbolic/primitives.py index 231ef872257cf888c1a5602e7548a7e6c34a761f..99d1e7f8eded6e20e8361d236ea152a4779b3950 100644 --- a/pytential/symbolic/primitives.py +++ b/pytential/symbolic/primitives.py @@ -23,7 +23,7 @@ THE SOFTWARE. """ import six -from six.moves import range, intern +from six.moves import intern from warnings import warn import numpy as np @@ -165,6 +165,7 @@ class Function(var): else: return var.__call__(self, operand) + real = Function("real") imag = Function("imag") conj = Function("conj") @@ -573,10 +574,12 @@ class IntGdSource(IntG): # {{{ geometric calculus -class _unspecified: + +class _unspecified: # noqa pass -def S(kernel, density, + +def S(kernel, density, # noqa qbx_forced_limit=_unspecified, source=None, target=None, kernel_arguments=None, **kwargs): @@ -641,178 +644,6 @@ def Dp(*args, **kwargs): # noqa # }}} -# {{{ differential operators on layer potentials - -def grad_S(kernel, arg, dim): - from pytools.obj_array import log_shape - arg_shape = log_shape(arg) - result = np.zeros(arg_shape+(dim,), dtype=object) - from pytools import indices_in_shape - for i in indices_in_shape(arg_shape): - for j in range(dim): - result[i+(j,)] = IntGdTarget(kernel, arg[i], j) - return result - - -def grad_D(kernel, arg, dim): - from pytools.obj_array import log_shape - arg_shape = log_shape(arg) - result = np.zeros(arg_shape+(dim,), dtype=object) - from pytools import indices_in_shape - for i in indices_in_shape(arg_shape): - for j in range(dim): - result[i+(j,)] = IntGdMixed(kernel, arg[i], j) - return result - - -def tangential_surf_grad_source_S(kernel, arg, dim=3): - from pytools.obj_array import make_obj_array - return make_obj_array([ - IntGdSource(kernel, arg, - ds_direction=make_tangent(i, dim, "src")) - for i in range(dim-1)]) - - -def surf_grad_S(kernel, arg, dim): - """ - :arg dim: The dimension of the ambient space. - """ - - return project_to_tangential(cse(grad_S(kernel, arg, dim))) - - -def div_S_volume(kernel, arg): - return sum(IntGdTarget(kernel, arg_n, n) for n, arg_n in enumerate(arg)) - - -def curl_S_volume(kernel, arg): - from pytools import levi_civita - from pytools.obj_array import make_obj_array - - return make_obj_array([ - sum( - levi_civita((l, m, n)) * IntGdTarget(kernel, arg[n], m) - for m in range(3) for n in range(3)) - for l in range(3)]) - - -def curl_curl_S_volume(k, arg): - # By vector identity, this is grad div S volume + k^2 S_k(arg), - # since S_k(arg) satisfies a Helmholtz equation. - - from pytools.obj_array import make_obj_array - - def swap_min_first(i, j): - if i < j: - return i, j - else: - return j, i - - return make_obj_array([ - sum(IntGd2Target(k, arg[m], *swap_min_first(m, n)) for m in range(3)) - for n in range(3)]) + k**2*S(k, arg) - - -def nxcurl_S(kernel, loc, arg): - """ - :arg loc: one of three values: - * +1 on the side of the surface toward which - the normal points ('exterior' of the surface), - * 0 on the surface, or evaluated at a volume target - * -1 on the interior of the surface. - """ - nxcurl_S = np.cross(normal(3), curl_S_volume(kernel, arg)) - assert loc in [-1, 0, 1], "invalid value for 'loc' (%s)" % loc - return nxcurl_S + loc*(1/2)*arg - - -def surface_laplacian_S_squared(u, invertibility_scale=0): - """ - :arg u: The field to which the surface Laplacian is applied. - """ - # http://wiki.tiker.net/HellsKitchen/SurfaceLaplacian - - Su = cse(S(0, u), "su_from_surflap") - - return ( - - 2*mean_curvature()*Sp(0, Su) - - ((Spp(0, Su)+Dp(0, Su))-(-1/4*u+Sp(0, Sp(0, u)))) - - invertibility_scale * mean(S(0, Su))*Ones()) - - -def S_surface_laplacian_S(u, dim, invertibility_scale=0, qbx_fix_scale=0): - """ - :arg u: The field to which the surface Laplacian is applied. - """ - - # This converges, but appears to work quite poorly compared to the above. - - tgrad_Su = cse( - project_to_tangential(grad_S(0, u, dim)), - "tgrad_su_from_surflap") - - return ( - - IntGdSource(0, Ones(), ds_direction=real(tgrad_Su)) - - 1j*IntGdSource(0, Ones(), ds_direction=imag(tgrad_Su)) - - invertibility_scale * S(0, Ones()*mean(S(0, u))) - - qbx_fix_scale * ( - u - # D+ - D- = identity (but QBX will compute the - # 'compact part of the identity' -- call that I*) - - ( - D(0, u, qbx_forced_limit=+1) - - D(0, u, qbx_forced_limit=-1)) - ) - # The above is I - I*, which means only the high-frequency - # bits of the identity are left. - ) - -# }}} - - -# {{{ geometric operations - -def xyz_to_tangential(xyz_vec, which=None): - d = len(xyz_vec) - x2l = xyz_to_local_matrix(d) - return np.dot(x2l[:-1], xyz_vec) - - -def tangential_to_xyz(tangential_vec, which=None): - d = len(tangential_vec) + 1 - x2l = xyz_to_local_matrix(d) - return np.dot(x2l[:-1].T, tangential_vec) - - -def project_to_tangential(xyz_vec, which=None): - return tangential_to_xyz( - cse(xyz_to_tangential(xyz_vec, which), which)) - - -def n_dot(vec, which=None): - return np.dot(normal(len(vec), which), vec) - - -def n_cross(vec, which=None): - nrm = normal(3, which) - - from pytools import levi_civita - from pytools.obj_array import make_obj_array - return make_obj_array([ - sum( - levi_civita((i, j, k)) * nrm[j] * vec[k] - for j in range(3) for k in range(3)) - for i in range(3)]) - - -def surf_n_cross(tangential_vec): - assert len(tangential_vec) == 2 - from pytools.obj_array import make_obj_array - return make_obj_array([-tangential_vec[1], tangential_vec[0]]) - -# }}} - - def pretty(expr): # Doesn't quite belong here, but this is exposed to the user as # "pytential.sym", so in here it goes. diff --git a/pytential/version.py b/pytential/version.py index 403f67f32928eda5924c824e5a4e652f86624059..d26fbc2f9341a880b1119e7e6079bd51e59e11b9 100644 --- a/pytential/version.py +++ b/pytential/version.py @@ -1,3 +1,2 @@ -VERSION = (2013, 1) +VERSION = (2016, 1) VERSION_TEXT = ".".join(str(i) for i in VERSION) - diff --git a/setup.cfg b/setup.cfg index 77f81d3bfbff0b437d7d65c0bcda8064887aa084..377c7dc2988b9234156f73afcaf748696205c754 100644 --- a/setup.cfg +++ b/setup.cfg @@ -1,3 +1,7 @@ [flake8] ignore = E126,E127,E128,E123,E226,E241,E242,E265,E402,W503 max-line-length=85 +exclude= + pytential/symbolic/pde/maxwell.py, + pytential/symbolic/old_diffop_primitives.py, + diff --git a/test/test_muller.py b/test/test_muller.py index 73de6b3df6dcb8e7dd5c0ed0501faa7c7f019943..30d1a47d8564f607a125bdfd0df5f6e67f811606 100644 --- a/test/test_muller.py +++ b/test/test_muller.py @@ -1,7 +1,4 @@ -from __future__ import division -from __future__ import absolute_import -from __future__ import print_function -from six.moves import range +from __future__ import division, absolute_import, print_function __copyright__ = "Copyright (C) 2014 Shidong Jiang, Andreas Kloeckner" @@ -25,6 +22,7 @@ OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE SOFTWARE. """ +from six.moves import range import numpy as np import pytest @@ -74,6 +72,7 @@ def fun1(z, n): return y + if __name__ == "__main__": import sys if len(sys.argv) > 1: diff --git a/test/test_tools.py b/test/test_tools.py index 7ca0d96c5290ccb30d9e415fd31106f63926eb53..0b003781c5b6fd28065193883bb891d60cfc18bb 100644 --- a/test/test_tools.py +++ b/test/test_tools.py @@ -29,14 +29,14 @@ import numpy.linalg as la def test_gmres(): n = 200 - A = ( + A = ( # noqa n * (np.eye(n) + 2j * np.eye(n)) + np.random.randn(n, n) + 1j * np.random.randn(n, n)) true_sol = np.random.randn(n) + 1j * np.random.randn(n) b = np.dot(A, true_sol) - A_func = lambda x: np.dot(A, x) + A_func = lambda x: np.dot(A, x) # noqa A_func.shape = A.shape A_func.dtype = A.dtype diff --git a/test/too_slow_test_helmholtz.py b/test/too_slow_test_helmholtz.py index 83f2f512a68bc5c9bc9b4044de6f1ddc63f36f0b..287a678007057a166daa0fe57b93a0b6a251da21 100644 --- a/test/too_slow_test_helmholtz.py +++ b/test/too_slow_test_helmholtz.py @@ -271,8 +271,8 @@ def run_dielectric_test(cl_ctx, queue, nelements, qbx_order, _, (H1_tgt_true,) = pot_p2p(queue, targets_1, h_sources_1, [h_strengths_1], out_host=True, k=K1) - err_F0_total = 0 - err_F1_total = 0 + err_F0_total = 0 # noqa + err_F1_total = 0 # noqa i_field = 0 @@ -287,11 +287,11 @@ def run_dielectric_test(cl_ctx, queue, nelements, qbx_order, continue if field_kind == pde_op.field_kind_e: - F0_tgt_true = E0_tgt_true - F1_tgt_true = E1_tgt_true + F0_tgt_true = E0_tgt_true # noqa + F1_tgt_true = E1_tgt_true # noqa elif field_kind == pde_op.field_kind_h: - F0_tgt_true = H0_tgt_true - F1_tgt_true = H1_tgt_true + F0_tgt_true = H0_tgt_true # noqa + F1_tgt_true = H1_tgt_true # noqa else: assert False @@ -301,8 +301,8 @@ def run_dielectric_test(cl_ctx, queue, nelements, qbx_order, rel_err_F0 = abs_err_F0/vec_norm(F0_tgt_true) # noqa rel_err_F1 = abs_err_F1/vec_norm(F1_tgt_true) # noqa - err_F0_total = max(rel_err_F0, err_F0_total) - err_F1_total = max(rel_err_F1, err_F1_total) + err_F0_total = max(rel_err_F0, err_F0_total) # noqa + err_F1_total = max(rel_err_F1, err_F1_total) # noqa print("Abs Err %s0" % field_kind_to_string(field_kind), abs_err_F0) print("Abs Err %s1" % field_kind_to_string(field_kind), abs_err_F1)