discretize a BilinearForm with spaces of different degree

psydac/api/discretization.py

@@ -29,7 +29,7 @@
 from psydac.api.fem          import DiscreteLinearForm
 from psydac.api.fem          import DiscreteFunctional
 from psydac.api.fem          import DiscreteSumForm
-from psydac.api.feec         import DiscreteDerham
+from psydac.api.feec         import DiscreteDerham, DiscreteDerhamMultipatch
 from psydac.api.glt          import DiscreteGltExpr
 from psydac.api.expr         import DiscreteExpr
 from psydac.api.equation     import DiscreteEquation
@@ -47,6 +47,7 @@
 __all__ = (
     'discretize',
     'discretize_derham',
+    'discretize_derham_multipatch',
     'reduce_space_degrees',
     'discretize_space',
     'discretize_domain'
@@ -175,11 +176,9 @@ def discretize_derham(derham, domain_h, *, get_H1vec_space=False, **kwargs):
     See Also
     --------
     discretize_space
-
     """
 
     ldim    = derham.shape
-    mapping = domain_h.domain.mapping # NOTE: assuming single-patch domain!
     bases   = ['B'] + ldim * ['M']
     spaces  = [discretize_space(V, domain_h, basis=basis, **kwargs)
                for V, basis in zip(derham.spaces, bases)]
@@ -192,7 +191,48 @@ def discretize_derham(derham, domain_h, *, get_H1vec_space=False, **kwargs):
         #We still need to specify the symbolic space because of "_recursive_element_of" not implemented in sympde
         spaces.append(Xh)
 
-    return DiscreteDerham(mapping, *spaces)
+    return DiscreteDerham(domain_h, *spaces)
+
+#==============================================================================
+def discretize_derham_multipatch(derham, domain_h, *args, **kwargs):
+    """
+    Create a discrete multipatch De Rham sequence from a symbolic one.
+
+    This function creates the broken discrete spaces from the symbolic ones, and then
+    creates a DiscreteDerhamMultipatch object from them.
+
+    Parameters
+    ----------
+    derham : sympde.topology.space.Derham
+        The symbolic Derham sequence.
+
+    domain_h : Geometry
+        Discrete domain where the spaces will be discretized.
+
+    **kwargs : dict
+        Optional parameters for the space discretization.
+
+    Returns
+    -------
+    DiscreteDerhamMultipatch
+      The discrete multipatch De Rham sequence containing the discrete spaces, 
+      differential operators and projectors.
+
+    See Also
+    --------
+    discretize_derham
+    discretize_space
+    """
+
+    ldim   = derham.shape
+    bases  = ['B'] + ldim * ['M']
+    spaces = [discretize_space(V, domain_h, *args, basis=basis, **kwargs) \
+            for V, basis in zip(derham.spaces, bases)]
+
+    return DiscreteDerhamMultipatch(
+        domain_h = domain_h,
+        spaces   = spaces
+    )
 
 #==============================================================================
 def reduce_space_degrees(V, Vh, *, basis='B', sequence='DR'):
@@ -329,7 +369,7 @@ def reduce_space_degrees(V, Vh, *, basis='B', sequence='DR'):
 
 #==============================================================================
 # TODO knots
-def discretize_space(V, domain_h, *, degree=None, multiplicity=None, knots=None, basis='B', sequence='DR'):
+def discretize_space(V, domain_h, *, degree=None, multiplicity=None, knots=None, basis='B', sequence='DR', pads=None):
     """
     This function creates the discretized space starting from the symbolic space.
 
@@ -457,6 +497,9 @@ def discretize_space(V, domain_h, *, degree=None, multiplicity=None, knots=None,
             min_coords = interior.min_coords
             max_coords = interior.max_coords
 
+            if pads == None: 
+                pads = [None] * len(degree_i)
+
             assert len(ncells) == len(periodic) == len(degree_i)  == len(multiplicity_i) == len(min_coords) == len(max_coords)
 
             if knots is None:
@@ -465,10 +508,10 @@ def discretize_space(V, domain_h, *, degree=None, multiplicity=None, knots=None,
                          for xmin, xmax, ne in zip(min_coords, max_coords, ncells)]
 
                 # Create 1D finite element spaces and precompute quadrature data
-                spaces[i] = [SplineSpace( p, multiplicity=m, grid=grid , periodic=P) for p,m,grid,P in zip(degree_i, multiplicity_i,grids, periodic)]
+                spaces[i] = [SplineSpace( p, multiplicity=m, grid=grid , periodic=P, pads=pd) for p,m,grid,P,pd in zip(degree_i, multiplicity_i,grids, periodic, pads)]
             else:
                  # Create 1D finite element spaces and precompute quadrature data
-                spaces[i] = [SplineSpace( p, knots=T , periodic=P) for p,T, P in zip(degree_i, knots[interior.name], periodic)]
+                spaces[i] = [SplineSpace( p, knots=T , periodic=P, pads=pd) for p,T,P,pd in zip(degree_i, knots[interior.name], periodic, pads)]
 
 
         carts    = create_cart(ddms, spaces)
@@ -607,9 +650,12 @@ def discretize(a, *args, **kwargs):
     elif isinstance(a, BasicFunctionSpace):
         return discretize_space(a, *args, **kwargs)
 
-    elif isinstance(a, Derham):
+    elif isinstance(a, Derham) and not a.V0.is_broken:
         return discretize_derham(a, *args, **kwargs)
 
+    elif isinstance(a, Derham) and a.V0.is_broken:
+        return discretize_derham_multipatch(a, *args, **kwargs)
+
     elif isinstance(a, Domain):
         return discretize_domain(a, *args, **kwargs)