Merge pull request #1459 from july-liuzhi/develop

update Develop
weihuayi · Jan 16, 2025 · ac5de6e · ac5de6e
2 parents 0596257 + 4585b7a
commit ac5de6e
Show file tree

Hide file tree

Showing 2 changed files with 33 additions and 8 deletions.
diff --git a/fealpy/fem/curlcurl_integrator.py b/fealpy/fem/curlcurl_integrator.py
@@ -32,6 +32,7 @@ def to_global_dof(self, space: _FS) -> TensorLike:
 
     def assembly(self, space, index: Index = _S, 
                              cellmeasure=None, out=None):
+        coef = self.coef
         mesh = space.mesh
         NC = mesh.number_of_cells()
         GD = space.mesh.geo_dimension()
@@ -44,10 +45,13 @@ def assembly(self, space, index: Index = _S,
         bcs, ws = qf.get_quadrature_points_and_weights()
 
         cphi = space.curl_basis(bcs) #(NQ, NC, ldof)
-        if GD==2:
-            cphi = cphi[..., None]
-        A = self.coef*bm.einsum("cqli, cqdi, c, q->cld", cphi, cphi, cm, ws)
-        if out is None:
-            return A
-        else:
-            out += A
+        val = process_coef_func(coef, bcs=bcs, mesh=mesh, etype='cell', index=index)
+        return bilinear_integral(cphi, cphi, ws, cm, val)
+
+        # if GD==2:
+        #     cphi = cphi[..., None]
+        # A = self.coef*bm.einsum("cqli, cqdi, c, q->cld", cphi, cphi, cm, ws)
+        # if out is None:
+        #     return A
+        # else:
+        #     out += A
diff --git a/fealpy/functionspace/first_nedelec_fe_space_3d.py b/fealpy/functionspace/first_nedelec_fe_space_3d.py
@@ -327,7 +327,7 @@ def face_basis(self, bcs, index=_S):
         eldof = self.dof.number_of_local_dofs("edge")
         ldof = 3*eldof + fldof
         gdof = self.dof.number_of_global_dofs()
-        glambda = mesh.grad_face_lambda()
+        glambda = mesh.grad_lambda(TD =2)
         ledge = bm.array([[1, 2],
                          [2, 0],
                          [0, 1]],device=self.device)
@@ -550,6 +550,27 @@ def curl_matrix(self):
         A = bm.to_numpy(A).reshape(-1)
         B = csr_matrix((A, (I, J)), shape=(gdof, gdof))
         return B
+
+    def face_mass_matrix(self, c = 1, index=_S):
+        """
+        @brief (n \times u, n \times v)_{Gamma_{robin}}
+        @param c 系数, 现在只考虑了 c 是常数的情况
+        """
+        p = self.p
+        mesh = self.mesh
+        qf = mesh.quadrature_formula(p+2,"face")
+        bcs, ws = qf.get_quadrature_points_and_weights()
+        face2dof = self.dof.face_to_dof()[index]
+        fm = self.mesh.entity_measure("face", index=index)
+        fphi = self.face_basis(bcs, index=index)
+
+        EMc = bm.einsum('fqlg, fqmg, q, f->flm', fphi, fphi, ws, fm)
+
+        gdof = self.dof.number_of_global_dofs()
+        I = bm.broadcast_to(face2dof[:,:, None], EMc.shape)
+        J = bm.broadcast_to(face2dof[:, None, :], EMc.shape)
+        EM = c*csr_matrix((EMc.flat, (I.flat, J.flat)), shape=(gdof, gdof))
+        return EM
 
     def source_vector(self, f):
         mesh = self.mesh