diff --git a/psydac/feec/multipatch/examples/__pycache__/__init__.cpython-310.pyc b/psydac/feec/multipatch/examples/__pycache__/__init__.cpython-310.pyc new file mode 100644 index 000000000..468b07c84 Binary files /dev/null and b/psydac/feec/multipatch/examples/__pycache__/__init__.cpython-310.pyc differ diff --git a/psydac/feec/multipatch/examples/__pycache__/h1_source_pbms_conga_2d.cpython-310.pyc b/psydac/feec/multipatch/examples/__pycache__/h1_source_pbms_conga_2d.cpython-310.pyc new file mode 100644 index 000000000..b24d0cabc Binary files /dev/null and b/psydac/feec/multipatch/examples/__pycache__/h1_source_pbms_conga_2d.cpython-310.pyc differ diff --git a/psydac/feec/multipatch/examples/__pycache__/hcurl_eigen_pbms_conga_2d.cpython-310.pyc b/psydac/feec/multipatch/examples/__pycache__/hcurl_eigen_pbms_conga_2d.cpython-310.pyc new file mode 100644 index 000000000..3f17954a8 Binary files /dev/null and b/psydac/feec/multipatch/examples/__pycache__/hcurl_eigen_pbms_conga_2d.cpython-310.pyc differ diff --git a/psydac/feec/multipatch/examples/__pycache__/hcurl_source_pbms_conga_2d.cpython-310.pyc b/psydac/feec/multipatch/examples/__pycache__/hcurl_source_pbms_conga_2d.cpython-310.pyc new file mode 100644 index 000000000..530a8b6e7 Binary files /dev/null and b/psydac/feec/multipatch/examples/__pycache__/hcurl_source_pbms_conga_2d.cpython-310.pyc differ diff --git a/psydac/feec/multipatch/examples/__pycache__/ppc_test_cases.cpython-310.pyc b/psydac/feec/multipatch/examples/__pycache__/ppc_test_cases.cpython-310.pyc new file mode 100644 index 000000000..f7849a90c Binary files /dev/null and b/psydac/feec/multipatch/examples/__pycache__/ppc_test_cases.cpython-310.pyc differ diff --git a/psydac/feec/multipatch/examples/h1_source_pbms_conga_2d.py b/psydac/feec/multipatch/examples/h1_source_pbms_conga_2d.py index cb03e4ff1..92b2451ba 100644 --- a/psydac/feec/multipatch/examples/h1_source_pbms_conga_2d.py +++ b/psydac/feec/multipatch/examples/h1_source_pbms_conga_2d.py @@ -25,6 +25,7 @@ from psydac.feec.multipatch.multipatch_domain_utilities import build_multipatch_domain from psydac.feec.multipatch.examples.ppc_test_cases import get_source_and_solution from psydac.feec.multipatch.utilities import time_count +from psydac.feec.multipatch.non_matching_operators import construct_V0_conforming_projection, construct_V1_conforming_projection from psydac.linalg.utilities import array_to_psydac from psydac.fem.basic import FemField @@ -92,7 +93,7 @@ def solve_h1_source_pbm( print('building the symbolic and discrete deRham sequences...') derham = Derham(domain, ["H1", "Hcurl", "L2"]) - derham_h = discretize(derham, domain_h, degree=degree, backend=PSYDAC_BACKENDS[backend_language]) + derham_h = discretize(derham, domain_h, degree=degree) # multi-patch (broken) spaces V0h = derham_h.V0 @@ -128,10 +129,8 @@ def solve_h1_source_pbm( print('conforming projection operators...') # conforming Projections (should take into account the boundary conditions of the continuous deRham sequence) - cP0 = derham_h.conforming_projection(space='V0', hom_bc=True, backend_language=backend_language) - cP0_m = cP0.to_sparse_matrix() - # cP1 = derham_h.conforming_projection(space='V1', hom_bc=True, backend_language=backend_language) - # cP1_m = cP1.to_sparse_matrix() + cP0_m = construct_V0_conforming_projection(V0h, domain_h, hom_bc=True) + # cP1_m = construct_V1_conforming_projection(V1h, domain_h, hom_bc=True) if not os.path.exists(plot_dir): os.makedirs(plot_dir) @@ -141,7 +140,7 @@ def lift_u_bc(u_bc): print('lifting the boundary condition in V0h... [warning: Not Tested Yet!]') # note: for simplicity we apply the full P1 on u_bc, but we only need to set the boundary dofs u_bc = lambdify(domain.coordinates, u_bc) - u_bc_log = [pull_2d_h1(u_bc, m) for m in mappings_list] + u_bc_log = [pull_2d_h1(u_bc, m.get_callable_mapping()) for m in mappings_list] # it's a bit weird to apply P1 on the list of (pulled back) logical fields -- why not just apply it on u_bc ? uh_bc = P0(u_bc_log) ubc_c = uh_bc.coeffs.toarray() @@ -176,7 +175,7 @@ def lift_u_bc(u_bc): if source_proj == 'P_geom': print('projecting the source with commuting projection P0...') f = lambdify(domain.coordinates, f_scal) - f_log = [pull_2d_h1(f, m) for m in mappings_list] + f_log = [pull_2d_h1(f, m.get_callable_mapping()) for m in mappings_list] f_h = P0(f_log) f_c = f_h.coeffs.toarray() b_c = dH0_m.dot(f_c) @@ -186,7 +185,7 @@ def lift_u_bc(u_bc): v = element_of(V0h.symbolic_space, name='v') expr = f_scal * v l = LinearForm(v, integral(domain, expr)) - lh = discretize(l, domain_h, V0h, backend=PSYDAC_BACKENDS[backend_language]) + lh = discretize(l, domain_h, V0h) b = lh.assemble() b_c = b.toarray() if plot_source: diff --git a/psydac/feec/multipatch/examples/hcurl_eigen_pbms_conga_2d.py b/psydac/feec/multipatch/examples/hcurl_eigen_pbms_conga_2d.py index 759efce0e..ce719ea66 100644 --- a/psydac/feec/multipatch/examples/hcurl_eigen_pbms_conga_2d.py +++ b/psydac/feec/multipatch/examples/hcurl_eigen_pbms_conga_2d.py @@ -18,6 +18,7 @@ from psydac.feec.multipatch.multipatch_domain_utilities import build_multipatch_domain from psydac.feec.multipatch.plotting_utilities import plot_field from psydac.feec.multipatch.utilities import time_count +from psydac.feec.multipatch.non_matching_operators import construct_V0_conforming_projection, construct_V1_conforming_projection def hcurl_solve_eigen_pbm(nc=4, deg=4, domain_name='pretzel_f', backend_language='python', mu=1, nu=1, gamma_h=10, sigma=None, nb_eigs=4, nb_eigs_plot=4, @@ -71,7 +72,7 @@ def hcurl_solve_eigen_pbm(nc=4, deg=4, domain_name='pretzel_f', backend_language print('building symbolic and discrete derham sequences...') derham = Derham(domain, ["H1", "Hcurl", "L2"]) - derham_h = discretize(derham, domain_h, degree=degree, backend=PSYDAC_BACKENDS[backend_language]) + derham_h = discretize(derham, domain_h, degree=degree) V0h = derham_h.V0 V1h = derham_h.V1 @@ -102,10 +103,8 @@ def hcurl_solve_eigen_pbm(nc=4, deg=4, domain_name='pretzel_f', backend_language print('conforming projection operators...') # conforming Projections (should take into account the boundary conditions of the continuous deRham sequence) - cP0 = derham_h.conforming_projection(space='V0', hom_bc=True, backend_language=backend_language, load_dir=m_load_dir) - cP1 = derham_h.conforming_projection(space='V1', hom_bc=True, backend_language=backend_language, load_dir=m_load_dir) - cP0_m = cP0.to_sparse_matrix() - cP1_m = cP1.to_sparse_matrix() + cP0_m = construct_V0_conforming_projection(V0h, domain_h, True) + cP1_m = construct_V1_conforming_projection(V1h, domain_h, True) print('broken differential operators...') bD0, bD1 = derham_h.broken_derivatives_as_operators @@ -214,10 +213,10 @@ def get_eigenvalues(nb_eigs, sigma, A_m, M_m): nc = 8 deg = 4 - domain_name = 'pretzel_f' - # domain_name = 'curved_L_shape' + #domain_name = 'pretzel_f' + domain_name = 'curved_L_shape' nc = 10 - deg = 2 + deg = 3 m_load_dir = 'matrices_{}_nc={}_deg={}/'.format(domain_name, nc, deg) run_dir = 'eigenpbm_{}_nc={}_deg={}/'.format(domain_name, nc, deg) diff --git a/psydac/feec/multipatch/examples/hcurl_source_pbms_conga_2d.py b/psydac/feec/multipatch/examples/hcurl_source_pbms_conga_2d.py index 8809898fe..55dd506e6 100644 --- a/psydac/feec/multipatch/examples/hcurl_source_pbms_conga_2d.py +++ b/psydac/feec/multipatch/examples/hcurl_source_pbms_conga_2d.py @@ -29,6 +29,8 @@ from psydac.linalg.utilities import array_to_psydac from psydac.fem.basic import FemField +from psydac.feec.multipatch.non_matching_operators import construct_V0_conforming_projection, construct_V1_conforming_projection + def solve_hcurl_source_pbm( nc=4, deg=4, domain_name='pretzel_f', backend_language=None, source_proj='P_geom', source_type='manu_J', eta=-10., mu=1., nu=1., gamma_h=10., @@ -107,7 +109,7 @@ def solve_hcurl_source_pbm( t_stamp = time_count(t_stamp) print('building discrete derham sequence...') - derham_h = discretize(derham, domain_h, degree=degree, backend=PSYDAC_BACKENDS[backend_language]) + derham_h = discretize(derham, domain_h, degree=degree) t_stamp = time_count(t_stamp) print('building commuting projection operators...') @@ -162,10 +164,8 @@ def solve_hcurl_source_pbm( t_stamp = time_count(t_stamp) print('building the conforming Projection operators and matrices...') # conforming Projections (should take into account the boundary conditions of the continuous deRham sequence) - cP0 = derham_h.conforming_projection(space='V0', hom_bc=True, backend_language=backend_language, load_dir=m_load_dir) - cP1 = derham_h.conforming_projection(space='V1', hom_bc=True, backend_language=backend_language, load_dir=m_load_dir) - cP0_m = cP0.to_sparse_matrix() - cP1_m = cP1.to_sparse_matrix() + cP0_m = construct_V0_conforming_projection(V0h, domain_h, hom_bc=True) + cP1_m = construct_V1_conforming_projection(V1h, domain_h, hom_bc=True) t_stamp = time_count(t_stamp) print('building the broken differential operators and matrices...') @@ -183,7 +183,7 @@ def lift_u_bc(u_bc): # note: for simplicity we apply the full P1 on u_bc, but we only need to set the boundary dofs u_bc_x = lambdify(domain.coordinates, u_bc[0]) u_bc_y = lambdify(domain.coordinates, u_bc[1]) - u_bc_log = [pull_2d_hcurl([u_bc_x, u_bc_y], m) for m in mappings_list] + u_bc_log = [pull_2d_hcurl([u_bc_x, u_bc_y], m.get_callable_mapping()) for m in mappings_list] # it's a bit weird to apply P1 on the list of (pulled back) logical fields -- why not just apply it on u_bc ? uh_bc = P1(u_bc_log) ubc_c = uh_bc.coeffs.toarray() @@ -240,7 +240,7 @@ def lift_u_bc(u_bc): print('projecting the source with commuting projection...') f_x = lambdify(domain.coordinates, f_vect[0]) f_y = lambdify(domain.coordinates, f_vect[1]) - f_log = [pull_2d_hcurl([f_x, f_y], m) for m in mappings_list] + f_log = [pull_2d_hcurl([f_x, f_y], m.get_callable_mapping()) for m in mappings_list] f_h = P1(f_log) f_c = f_h.coeffs.toarray() b_c = dH1_m.dot(f_c) @@ -251,7 +251,7 @@ def lift_u_bc(u_bc): v = element_of(V1h.symbolic_space, name='v') expr = dot(f_vect,v) l = LinearForm(v, integral(domain, expr)) - lh = discretize(l, domain_h, V1h, backend=PSYDAC_BACKENDS[backend_language]) + lh = discretize(l, domain_h, V1h) b = lh.assemble() b_c = b.toarray() if plot_source: diff --git a/psydac/feec/multipatch/examples/mixed_source_pbms_conga_2d.py b/psydac/feec/multipatch/examples/mixed_source_pbms_conga_2d.py index 5fba55dc3..81212af59 100644 --- a/psydac/feec/multipatch/examples/mixed_source_pbms_conga_2d.py +++ b/psydac/feec/multipatch/examples/mixed_source_pbms_conga_2d.py @@ -28,6 +28,8 @@ from psydac.feec.multipatch.examples.hcurl_eigen_pbms_conga_2d import get_eigenvalues from psydac.feec.multipatch.utilities import time_count +from psydac.feec.multipatch.non_matching_operators import construct_V0_conforming_projection, construct_V1_conforming_projection + def solve_magnetostatic_pbm( nc=4, deg=4, domain_name='pretzel_f', backend_language=None, source_proj='P_L2_wcurl_J', source_type='dipole_J', bc_type='metallic', @@ -120,7 +122,7 @@ def solve_magnetostatic_pbm( print('building symbolic and discrete derham sequences...') derham = Derham(domain, ["H1", "Hcurl", "L2"]) - derham_h = discretize(derham, domain_h, degree=degree, backend=PSYDAC_BACKENDS[backend_language]) + derham_h = discretize(derham, domain_h, degree=degree) V0h = derham_h.V0 V1h = derham_h.V1 @@ -137,18 +139,18 @@ def solve_magnetostatic_pbm( # these physical projection operators should probably be in the interface... def P0_phys(f_phys): f = lambdify(domain.coordinates, f_phys) - f_log = [pull_2d_h1(f, m) for m in mappings_list] + f_log = [pull_2d_h1(f, m.get_callable_mapping()) for m in mappings_list] return P0(f_log) def P1_phys(f_phys): f_x = lambdify(domain.coordinates, f_phys[0]) f_y = lambdify(domain.coordinates, f_phys[1]) - f_log = [pull_2d_hcurl([f_x, f_y], m) for m in mappings_list] + f_log = [pull_2d_hcurl([f_x, f_y], m.get_callable_mapping()) for m in mappings_list] return P1(f_log) def P2_phys(f_phys): f = lambdify(domain.coordinates, f_phys) - f_log = [pull_2d_l2(f, m) for m in mappings_list] + f_log = [pull_2d_l2(f, m.get_callable_mapping()) for m in mappings_list] return P2(f_log) I0_m = IdLinearOperator(V0h).to_sparse_matrix() @@ -175,10 +177,9 @@ def P2_phys(f_phys): print('conforming projection operators...') # conforming Projections (should take into account the boundary conditions of the continuous deRham sequence) - cP0 = derham_h.conforming_projection(space='V0', hom_bc=hom_bc, backend_language=backend_language, load_dir=m_load_dir) - cP1 = derham_h.conforming_projection(space='V1', hom_bc=hom_bc, backend_language=backend_language, load_dir=m_load_dir) - cP0_m = cP0.to_sparse_matrix() - cP1_m = cP1.to_sparse_matrix() + cP0_m = construct_V0_conforming_projection(V0h, domain_h, hom_bc=True) + cP1_m = construct_V1_conforming_projection(V1h, domain_h, hom_bc=True) + print('broken differential operators...') bD0, bD1 = derham_h.broken_derivatives_as_operators diff --git a/psydac/feec/multipatch/tests/__init__.py b/psydac/feec/multipatch/tests/__init__.py new file mode 100644 index 000000000..e69de29bb diff --git a/psydac/feec/multipatch/tests/test_feec_poisson_multipatch_2d.py b/psydac/feec/multipatch/tests/test_feec_poisson_multipatch_2d.py index 8eea88ec0..7fd6eb040 100644 --- a/psydac/feec/multipatch/tests/test_feec_poisson_multipatch_2d.py +++ b/psydac/feec/multipatch/tests/test_feec_poisson_multipatch_2d.py @@ -16,7 +16,7 @@ def test_poisson_pretzel_f(): backend_language='pyccel-gcc', plot_source=False, plot_dir='./plots/h1_tests_source_february/'+run_dir) - + print(l2_error) assert abs(l2_error-8.054935880021907e-05)<1e-10 #============================================================================== @@ -30,3 +30,6 @@ def teardown_module(): def teardown_function(): from sympy.core import cache cache.clear_cache() + +if __name__ == '__main__': + test_poisson_pretzel_f()