Skip to content

Commit a3ac0f9

Browse files
jorgensdjorgensd
authored
Address changes from: FEniCS/dolfinx#3740 (#184)
Requires: FEniCS/dolfinx#3877
1 parent f492bcc commit a3ac0f9

6 files changed

Lines changed: 46 additions & 51 deletions

File tree

‎cpp/assemble_matrix.cpp‎

Lines changed: 5 additions & 3 deletions
Original file line numberDiff line numberDiff line change
@@ -610,7 +610,7 @@ void assemble_matrix_impl(
610610
cell_info0 = std::span(mesh0->topology()->get_cell_permutation_info());
611611
cell_info1 = std::span(mesh1->topology()->get_cell_permutation_info());
612612
}
613-
for (int i : a.integral_ids(dolfinx::fem::IntegralType::cell))
613+
for (int i = 0; i < a.num_integrals(dolfinx::fem::IntegralType::cell, 0); ++i)
614614
{
615615
const auto& fn = a.kernel(dolfinx::fem::IntegralType::cell, i, 0);
616616
const auto& [coeffs, cstride]
@@ -628,7 +628,9 @@ void assemble_matrix_impl(
628628
cstride, constants, cell_info0, cell_info1, mpc0, mpc1);
629629
}
630630

631-
for (int i : a.integral_ids(dolfinx::fem::IntegralType::exterior_facet))
631+
for (int i = 0;
632+
i < a.num_integrals(dolfinx::fem::IntegralType::exterior_facet, 0); ++i)
633+
632634
{
633635
const auto& fn = a.kernel(dolfinx::fem::IntegralType::exterior_facet, i, 0);
634636
const auto& [coeffs, cstride]
@@ -646,7 +648,7 @@ void assemble_matrix_impl(
646648
cstride, constants, cell_info0, cell_info1, mpc0, mpc1);
647649
}
648650

649-
if (a.integral_ids(dolfinx::fem::IntegralType::interior_facet).size() > 0)
651+
if (a.num_integrals(dolfinx::fem::IntegralType::interior_facet, 0) > 0)
650652
throw std::runtime_error("Not implemented yet");
651653
}
652654
//-----------------------------------------------------------------------------

‎cpp/assemble_vector.cpp‎

Lines changed: 4 additions & 3 deletions
Original file line numberDiff line numberDiff line change
@@ -145,7 +145,7 @@ void _assemble_vector(
145145

146146
const auto fetch_cell
147147
= [&](std::span<const std::int32_t> entity) { return entity.front(); };
148-
for (int i : L.integral_ids(dolfinx::fem::IntegralType::cell))
148+
for (int i = 0; i < L.num_integrals(dolfinx::fem::IntegralType::cell, 0); ++i)
149149
{
150150
const auto& coeffs = coefficients.at({dolfinx::fem::IntegralType::cell, i});
151151

@@ -189,7 +189,8 @@ void _assemble_vector(
189189
// Prepare permutations for exterior and interior facet integrals
190190

191191
// Assemble exterior facet integral kernels
192-
for (int i : L.integral_ids(dolfinx::fem::IntegralType::exterior_facet))
192+
for (int i = 0;
193+
i < L.num_integrals(dolfinx::fem::IntegralType::exterior_facet, 0); ++i)
193194
{
194195
const auto& fn = L.kernel(dolfinx::fem::IntegralType::exterior_facet, i, 0);
195196
const auto& coeffs
@@ -234,7 +235,7 @@ void _assemble_vector(
234235
assemble_local_exterior_facet_vector);
235236
}
236237

237-
if (L.integral_ids(dolfinx::fem::IntegralType::interior_facet).size() > 0)
238+
if (L.num_integrals(dolfinx::fem::IntegralType::interior_facet, 0) > 0)
238239
{
239240
throw std::runtime_error(
240241
"Interior facet integrals currently not supported");

‎cpp/lifting.h‎

Lines changed: 13 additions & 11 deletions
Original file line numberDiff line numberDiff line change
@@ -235,7 +235,8 @@ void apply_lifting(
235235

236236
const auto fetch_cells
237237
= [&](std::span<const std::int32_t> entity) { return entity.front(); };
238-
for (int i : a->integral_ids(dolfinx::fem::IntegralType::cell))
238+
for (int i = 0; i < a->num_integrals(dolfinx::fem::IntegralType::cell, 0);
239+
++i)
239240
{
240241
const auto& coeffs = coefficients.at({dolfinx::fem::IntegralType::cell, i});
241242
const auto& kernel = a->kernel(dolfinx::fem::IntegralType::cell, i, 0);
@@ -307,7 +308,8 @@ void apply_lifting(
307308
}
308309

309310
// Get number of cells per facet to be able to get the facet permutation
310-
for (int i : a->integral_ids(dolfinx::fem::IntegralType::exterior_facet))
311+
for (int i = 0;
312+
i < a->num_integrals(dolfinx::fem::IntegralType::exterior_facet, 0); ++i)
311313
{
312314
const auto& coeffs
313315
= coefficients.at({dolfinx::fem::IntegralType::exterior_facet, i});
@@ -388,7 +390,7 @@ void apply_lifting(
388390
b, active_facets, active_facets0, active_facets1, *dofmap0, *dofmap1,
389391
bc_values1, bc_markers1, *mpc0, fetch_cells, lift_bc_exterior_facet);
390392
}
391-
if (a->integral_ids(dolfinx::fem::IntegralType::interior_facet).size() > 0)
393+
if (a->num_integrals(dolfinx::fem::IntegralType::interior_facet, 0) > 0)
392394
{
393395

394396
throw std::runtime_error(
@@ -518,25 +520,24 @@ void apply_lifting(
518520
"Mismatch in size between a and bcs in assembler.");
519521
}
520522
for (std::size_t j = 0; j < a.size(); ++j)
521-
{
523+
{
522524
if (a[j] and !bcs1[j].empty())
523525
{
524-
if (x0.empty())
526+
if (x0.empty())
525527
{
526528
impl::apply_lifting<std::complex<double>>(
527529
b, a[j], bcs1[j], std::span<const std::complex<double>>(), scale,
528530
mpc);
529531
}
530532
else
531533
{
532-
impl::apply_lifting<std::complex<double>>(b, a[j], bcs1[j], x0[j], scale,
533-
mpc);
534+
impl::apply_lifting<std::complex<double>>(b, a[j], bcs1[j], x0[j],
535+
scale, mpc);
534536
}
535537
}
536538
}
537539
}
538540

539-
540541
/// Modify b such that:
541542
///
542543
/// b <- b - scale * K^T (A_j (g_j 0 x0_j))
@@ -584,7 +585,7 @@ void apply_lifting(
584585
if (x0.empty())
585586
{
586587
impl::apply_lifting<float>(b, a[j], bcs1[j], std::span<const float>(),
587-
scale, mpc);
588+
scale, mpc);
588589
}
589590
else
590591
{
@@ -645,12 +646,13 @@ void apply_lifting(
645646
if (x0.empty())
646647
{
647648
impl::apply_lifting<std::complex<float>>(
648-
b, a[j], bcs1[j], std::span<const std::complex<float>>(), scale, mpc);
649+
b, a[j], bcs1[j], std::span<const std::complex<float>>(), scale,
650+
mpc);
649651
}
650652
else
651653
{
652654
impl::apply_lifting<std::complex<float>>(b, a[j], bcs1[j], x0[j], scale,
653-
mpc);
655+
mpc);
654656
}
655657
}
656658
}

‎cpp/utils.h‎

Lines changed: 8 additions & 10 deletions
Original file line numberDiff line numberDiff line change
@@ -313,35 +313,33 @@ void build_standard_pattern(dolfinx::la::SparsityPattern& pattern,
313313

314314
for (auto type : types)
315315
{
316-
std::vector<int> ids = a.integral_ids(type);
317316
switch (type)
318317
{
319318
case dolfinx::fem::IntegralType::cell:
320-
for (int id : ids)
319+
for (int i = 0; i < a.num_integrals(type, 0); ++i)
321320
{
322321
dolfinx::fem::sparsitybuild::cells(
323-
pattern,
324-
{a.domain_arg(type, 0, id, 0), a.domain_arg(type, 1, id, 0)},
322+
pattern, {a.domain_arg(type, 0, i, 0), a.domain_arg(type, 1, i, 0)},
325323
{{dofmaps[0], dofmaps[1]}});
326324
}
327325
break;
328326
case dolfinx::fem::IntegralType::interior_facet:
329-
for (int id : ids)
327+
for (int i = 0; i < a.num_integrals(type, 0); ++i)
330328
{
331329
dolfinx::fem::sparsitybuild::interior_facets(
332330
pattern,
333-
{extract_cells(a.domain_arg(type, 0, id, 0)),
334-
extract_cells(a.domain_arg(type, 1, id, 0))},
331+
{extract_cells(a.domain_arg(type, 0, i, 0)),
332+
extract_cells(a.domain_arg(type, 1, i, 0))},
335333
{{dofmaps[0], dofmaps[1]}});
336334
}
337335
break;
338336
case dolfinx::fem::IntegralType::exterior_facet:
339-
for (int id : ids)
337+
for (int i = 0; i < a.num_integrals(type, 0); ++i)
340338
{
341339
dolfinx::fem::sparsitybuild::cells(
342340
pattern,
343-
{extract_cells(a.domain_arg(type, 0, id, 0)),
344-
extract_cells(a.domain_arg(type, 1, id, 0))},
341+
{extract_cells(a.domain_arg(type, 0, i, 0)),
342+
extract_cells(a.domain_arg(type, 1, i, 0))},
345343
{{dofmaps[0], dofmaps[1]}});
346344
}
347345
break;

‎python/src/dolfinx_mpc/numba/assemble_matrix.py‎

Lines changed: 8 additions & 11 deletions
Original file line numberDiff line numberDiff line change
@@ -104,8 +104,7 @@ def assemble_matrix(
104104
tdim = V.mesh.topology.dim
105105

106106
# Assemble over cells
107-
subdomain_ids = form._cpp_object.integral_ids(_fem.IntegralType.cell)
108-
num_cell_integrals = len(subdomain_ids)
107+
num_cell_integrals = form._cpp_object.num_integrals(_fem.IntegralType.cell, 0)
109108

110109
e0 = form.function_spaces[0].element
111110
e1 = form.function_spaces[1].element
@@ -139,11 +138,10 @@ def assemble_matrix(
139138
# NOTE: This depends on enum ordering in ufcx.h
140139
cell_form_pos = ufcx_form.form_integral_offsets[0]
141140
V.mesh.topology.create_entity_permutations()
142-
for i, id in enumerate(subdomain_ids):
143-
coeffs_i = form_coeffs[(_fem.IntegralType.cell, id)]
144-
141+
for i in range(num_cell_integrals):
142+
coeffs_i = form_coeffs[(_fem.IntegralType.cell, i)]
145143
cell_kernel = getattr(ufcx_form.form_integrals[cell_form_pos + i], f"tabulate_tensor_{nptype}")
146-
active_cells = form._cpp_object.domains(_fem.IntegralType.cell, id)
144+
active_cells = form._cpp_object.domains(_fem.IntegralType.cell, i)
147145
assemble_slave_cells(
148146
A.handle,
149147
cell_kernel,
@@ -160,8 +158,7 @@ def assemble_matrix(
160158
)
161159

162160
# Assemble over exterior facets
163-
subdomain_ids = form._cpp_object.integral_ids(_fem.IntegralType.exterior_facet)
164-
num_exterior_integrals = len(subdomain_ids)
161+
num_exterior_integrals = form.num_integrals(_fem.IntegralType.exterior_facet, 0)
165162

166163
if num_exterior_integrals > 0:
167164
V.mesh.topology.create_entities(tdim - 1)
@@ -175,10 +172,10 @@ def assemble_matrix(
175172
perm = (cell_perms, form._cpp_object.needs_facet_permutations, facet_perms)
176173
# NOTE: This depends on enum ordering in ufcx.h
177174
ext_facet_pos = ufcx_form.form_integral_offsets[1]
178-
for i, id in enumerate(subdomain_ids):
175+
for i in range(num_exterior_integrals):
179176
facet_kernel = getattr(ufcx_form.form_integrals[ext_facet_pos + i], f"tabulate_tensor_{nptype}")
180-
facets = form._cpp_object.domains(_fem.IntegralType.exterior_facet, id)
181-
coeffs_i = form_coeffs[(_fem.IntegralType.exterior_facet, id)]
177+
facets = form._cpp_object.domains(_fem.IntegralType.exterior_facet, i)
178+
coeffs_i = form_coeffs[(_fem.IntegralType.exterior_facet, i)]
182179
facet_info = pack_slave_facet_info(facets, slave_cells)
183180
num_facets_per_cell = len(V.mesh.topology.connectivity(tdim, tdim - 1).links(0))
184181
assemble_exterior_slave_facets(

‎python/src/dolfinx_mpc/numba/assemble_vector.py‎

Lines changed: 8 additions & 13 deletions
Original file line numberDiff line numberDiff line change
@@ -95,9 +95,6 @@ def assemble_vector(form: _forms, constraint: MultiPointConstraint, b: Optional[
9595
if e0.needs_dof_transformations:
9696
raise NotImplementedError("Dof transformations not implemented")
9797
# Assemble over cells
98-
subdomain_ids = form._cpp_object.integral_ids(_fem.IntegralType.cell)
99-
num_cell_integrals = len(subdomain_ids)
100-
10198
if _PETSc.ScalarType == numpy.float32: # type: ignore
10299
nptype = "float32"
103100
elif _PETSc.ScalarType == numpy.float64: # type: ignore
@@ -109,15 +106,15 @@ def assemble_vector(form: _forms, constraint: MultiPointConstraint, b: Optional[
109106
else:
110107
raise RuntimeError(f"Unsupported scalar type {_PETSc.ScalarType}.") # type: ignore
111108
ufcx_form = form.ufcx_form
112-
if num_cell_integrals > 0:
109+
if (num_cell_integrals := form.num_integrals(_fem.IntegralType.cell, 0)) > 0:
113110
V.mesh.topology.create_entity_permutations()
114111

115112
# NOTE: This depends on enum ordering in ufcx.h
116113
cell_form_pos = ufcx_form.form_integral_offsets[0]
117-
for i, id in enumerate(subdomain_ids):
114+
for i in range(num_cell_integrals):
118115
cell_kernel = getattr(ufcx_form.form_integrals[cell_form_pos + i], f"tabulate_tensor_{nptype}")
119-
active_cells = form._cpp_object.domains(_fem.IntegralType.cell, id)
120-
coeffs_i = form_coeffs[(_fem.IntegralType.cell, id)]
116+
active_cells = form._cpp_object.domains(_fem.IntegralType.cell, i)
117+
coeffs_i = form_coeffs[(_fem.IntegralType.cell, i)]
121118
with vector.localForm() as b:
122119
assemble_cells(
123120
numpy.asarray(b),
@@ -134,9 +131,7 @@ def assemble_vector(form: _forms, constraint: MultiPointConstraint, b: Optional[
134131
)
135132

136133
# Assemble exterior facet integrals
137-
subdomain_ids = form._cpp_object.integral_ids(_fem.IntegralType.exterior_facet)
138-
num_exterior_integrals = len(subdomain_ids)
139-
if num_exterior_integrals > 0:
134+
if (num_exterior_integrals := form.num_integrals(_fem.IntegralType.exterior_facet, 0)) > 0:
140135
V.mesh.topology.create_entities(tdim - 1)
141136
V.mesh.topology.create_connectivity(tdim - 1, tdim)
142137
# Get facet permutations if required
@@ -146,10 +141,10 @@ def assemble_vector(form: _forms, constraint: MultiPointConstraint, b: Optional[
146141
perm = (cell_perms, form._cpp_object.needs_facet_permutations, facet_perms)
147142
# NOTE: This depends on enum ordering in ufcx.h
148143
ext_facet_pos = ufcx_form.form_integral_offsets[1]
149-
for i, id in enumerate(subdomain_ids):
144+
for i in range(num_exterior_integrals):
150145
facet_kernel = getattr(ufcx_form.form_integrals[ext_facet_pos + i], f"tabulate_tensor_{nptype}")
151-
coeffs_i = form_coeffs[(_fem.IntegralType.exterior_facet, id)]
152-
facets = form._cpp_object.domains(_fem.IntegralType.exterior_facet, id)
146+
coeffs_i = form_coeffs[(_fem.IntegralType.exterior_facet, i)]
147+
facets = form._cpp_object.domains(_fem.IntegralType.exterior_facet, i)
153148
facet_info = pack_slave_facet_info(facets, slave_cells)
154149
num_facets_per_cell = len(V.mesh.topology.connectivity(tdim, tdim - 1).links(0))
155150
with vector.localForm() as b:

0 commit comments

Comments
 (0)