diff --git a/hog/forms.py b/hog/forms.py index e8dc31b1c7c1926cfa900479aa5f3b0d52111645..8ae046393a7a8ad2b0330b12c960bb6eccfa7c1b 100644 --- a/hog/forms.py +++ b/hog/forms.py @@ -767,7 +767,6 @@ The resulting matrix must be multiplied with a vector of ones to be used as the docstring=docstring, ) - def divdiv( trial: TrialSpace, test: TestSpace, @@ -813,6 +812,43 @@ Weak formulation ) +def k_divdiv( + trial: TrialSpace, + test: TestSpace, + geometry: ElementGeometry, + symbolizer: Symbolizer, + coefficient_function_space: Optional[FunctionSpace] = None, + blending: GeometryMap = IdentityMap(), +) -> Form: + docstring = f""" +divdiv operator which is the compressible part of full Stokes operator + +Geometry map: {blending} + +Weak formulation + + u: trial function (vectorial space: {trial}) + v: test function (vectorial space: {test}) + + ∫ μ ( ∇ · u ) · ( ∇ · v ) +""" + + from hog.recipes.integrands.volume.divdiv import integrand + + return process_integrand( + integrand, + trial, + test, + geometry, + symbolizer, + blending=blending, + is_symmetric=trial == test, + docstring=docstring, + fe_coefficients={ + "k": coefficient_function_space + } + ) + def advection( trial: TrialSpace, test: TestSpace, diff --git a/hog/recipes/integrands/volume/divdiv.py b/hog/recipes/integrands/volume/divdiv.py index 80315d110b8eb18e74902c1f20009b27c6ca4f72..5a417a617fc132e47dfac2f416611bd6d2cf3c0f 100644 --- a/hog/recipes/integrands/volume/divdiv.py +++ b/hog/recipes/integrands/volume/divdiv.py @@ -26,8 +26,9 @@ def integrand( grad_u, grad_v, tabulate, + k, **_, ): div_u = (jac_b_inv.T * tabulate(jac_a_inv.T * grad_u)).trace() div_v = (jac_b_inv.T * tabulate(jac_a_inv.T * grad_v)).trace() - return div_u * div_v * tabulate(jac_a_abs_det) * jac_b_abs_det + return (k["k"] if "k" in k else 1.0) * div_u * div_v * tabulate(jac_a_abs_det) * jac_b_abs_det