diff --git a/example/experiment/fem/maxwell_eq_n1fem_dirichlet.py b/example/experiment/fem/maxwell_eq_n1fem_dirichlet.py index e1cd080d7..10df8b7d9 100644 --- a/example/experiment/fem/maxwell_eq_n1fem_dirichlet.py +++ b/example/experiment/fem/maxwell_eq_n1fem_dirichlet.py @@ -64,7 +64,9 @@ args = parser.parse_args() p = args.degree maxit = args.maxit +maxit = 3 dim = args.dim +dim = 3 backend = args.backend bm.set_backend(backend) if dim == 2: diff --git a/example/experiment/fem/maxwell_eq_n2fem_dirichlet.py b/example/experiment/fem/maxwell_eq_n2fem_dirichlet.py index f6fffe86b..3719ac329 100644 --- a/example/experiment/fem/maxwell_eq_n2fem_dirichlet.py +++ b/example/experiment/fem/maxwell_eq_n2fem_dirichlet.py @@ -65,9 +65,11 @@ args = parser.parse_args() p = args.degree maxit = args.maxit +maxit = 3 dim = args.dim +dim = 3 backend = args.backend -bm.set_backend("pytorch") +bm.set_backend(backend) if dim == 2: pde = PDE2d() else: @@ -79,8 +81,8 @@ '$||\\nabla \\times u - \\nabla_h \\times u_h||_0$ with k=2', '$||\\nabla \\times u - \\nabla_h \\times u_h||_0$ with k=3', '$||\\nabla \\times u - \\nabla_h \\times u_h||_0$ with k=4'] -errorMatrix = bm.zeros((len(errorType), maxit), dtype=bm.float64) -NDof = bm.zeros(maxit, dtype=bm.float64) +errorMatrix = np.zeros((len(errorType), maxit), dtype=np.float64) +NDof = np.zeros(maxit, dtype=np.float64) tmr = timer() next(tmr) diff --git a/fealpy/functionspace/first_nedelec_fe_space_3d.py b/fealpy/functionspace/first_nedelec_fe_space_3d.py index d0d2904f2..b1802915f 100644 --- a/fealpy/functionspace/first_nedelec_fe_space_3d.py +++ b/fealpy/functionspace/first_nedelec_fe_space_3d.py @@ -24,6 +24,7 @@ def __init__(self, mesh, p): self.multiindex = mesh.multi_index_matrix(p,3) self.ftype = mesh.ftype self.itype = mesh.itype + self.device = mesh.device def number_of_local_dofs(self, doftype='all'): p = self.p @@ -74,7 +75,7 @@ def face_to_dof(self): edge = self.mesh.entity('edge') face = self.mesh.entity('face') - f2d = bm.zeros((NF, fdof), dtype=self.itype) + f2d = bm.zeros((NF, fdof), device=self.device,dtype=self.itype) #f2d[:, :eldof*3] = e2dof[f2e].reshape(NF, eldof*3) f2d = bm.set_at(f2d,(slice(None),slice(eldof*3)),e2dof[f2e].reshape(NF, eldof*3)) s = [1, 0, 0] @@ -224,7 +225,7 @@ def basis(self, bcs, index=_S): c2esign = mesh.cell_to_edge_sign() #(NC, 3, 2) - l = bm.zeros((4, )+bcs[None,: ,0,None, None].shape, dtype=self.ftype) + l = bm.zeros((4, )+bcs[None,: ,0,None, None].shape,device=self.device, dtype=self.ftype) l = bm.set_at(l,(0),bcs[None, :,0,None, None]) l = bm.set_at(l,(1),bcs[None, :,1,None, None]) @@ -300,7 +301,7 @@ def face_internal_basis(self, bcs, index=_S): fdof = self.dof.number_of_local_dofs("face") glambda = mesh.grad_face_lambda() - val = bm.zeros((NF,) + bcs.shape[:-1]+(fdof, 3), dtype=self.ftype) + val = bm.zeros((NF,) + bcs.shape[:-1]+(fdof, 3),device=self.device, dtype=self.ftype) l = bm.zeros((3, )+bcs[None, :,0, None, None].shape, dtype=self.ftype) l = bm.set_at(l,(0),bcs[None, :,0, None, None]) @@ -338,7 +339,7 @@ def face_basis(self, bcs, index=_S): # edge basis phi = self.bspace.basis(bcs, p=p) multiIndex = self.mesh.multi_index_matrix(p, 2) - val = bm.zeros((NF,) + bcs.shape[:-1]+(ldof, 3), dtype=self.ftype) + val = bm.zeros((NF,) + bcs.shape[:-1]+(ldof, 3),device=self.device, dtype=self.ftype) for i in range(3): phie = phi[:, :, multiIndex[:, i]==0] c2esi = c2esign[:, i] @@ -384,7 +385,7 @@ def curl_basis(self, bcs): phi = self.bspace.basis(bcs, p=p) gphi = self.bspace.grad_basis(bcs, p=p) multiIndex = self.mesh.multi_index_matrix(p, 3) - val = bm.zeros((NC,)+bcs.shape[:-1]+(ldof, 3), dtype=self.ftype) + val = bm.zeros((NC,)+bcs.shape[:-1]+(ldof, 3),device=self.device, dtype=self.ftype) # edge basis locEdgeDual = bm.tensor([[2, 3], [1, 3], [1, 2], [0, 3], [0, 2], [0, 1]]) diff --git a/fealpy/functionspace/second_nedelec_fe_space_2d.py b/fealpy/functionspace/second_nedelec_fe_space_2d.py index c66edd5d6..8da8e07cf 100644 --- a/fealpy/functionspace/second_nedelec_fe_space_2d.py +++ b/fealpy/functionspace/second_nedelec_fe_space_2d.py @@ -96,7 +96,7 @@ def cell_to_dof(self): c2e = self.mesh.cell_to_edge() c2esign = self.mesh.cell_to_face_sign() - c2d = bm.zeros((NC, ldof), dtype=self.itype) + c2d = bm.zeros((NC, ldof), dtype=bm.int64) # c2d[:, e2ldof] : (NC, 3, p+1), e2dof[c2e] : (NC, 3, p+1) tmp = e2dof[c2e] # tmp[~c2esign] = tmp[~c2esign, ::-1] diff --git a/fealpy/functionspace/second_nedelec_fe_space_3d.py b/fealpy/functionspace/second_nedelec_fe_space_3d.py index c54402f6b..883ffea9e 100644 --- a/fealpy/functionspace/second_nedelec_fe_space_3d.py +++ b/fealpy/functionspace/second_nedelec_fe_space_3d.py @@ -26,7 +26,7 @@ def __init__(self, mesh, p): self.multiindex3 = mesh.multi_index_matrix(p,3) self.ftype = mesh.ftype self.itype = mesh.itype - + self.device=mesh.device def edge_to_local_face_dof(self): multiindex = self.multiindex2 @@ -53,7 +53,7 @@ def face_to_local_dof(self): ldof = self.number_of_local_dofs() fdof = self.number_of_local_dofs('faceall') - f2ld = bm.zeros((4, fdof), dtype=self.itype) + f2ld = bm.zeros((4, fdof),device=self.device, dtype=self.itype) eldof = self.edge_to_local_face_dof() nldof = bm.tensor([[eldof[(i+1)%3, 0], eldof[(i+2)%3, -1]] for i in range(3)]) @@ -91,10 +91,10 @@ def cell_to_dof(self): cdof = self.number_of_local_dofs('cell') fdof = self.number_of_local_dofs('faceall') - isndof = bm.zeros(ldof, dtype=bm.bool) + isndof = bm.zeros(ldof, device=self.device, dtype=bm.bool) isndof[f2ld] = True - c2d = bm.zeros((NC, ldof), dtype=bm.int64)# + c2d = bm.zeros((NC, ldof), device=self.device, dtype=bm.int64)# idx = bm.zeros((NC, 3), dtype=self.itype) fe = bm.tensor([[0, 1], [0, 2], [1, 2]], dtype=self.itype) #局部边 for i in range(4): @@ -235,7 +235,7 @@ def basis(self, bc,index=_S): c2v = self.basis_vector() #(NC, ldof, GD) shape = bc.shape[:-1] - val = bm.zeros((NC,)+ shape + (ldof, GD), dtype=self.ftype) + val = bm.zeros((NC,)+ shape + (ldof, GD),device=self.device, dtype=self.ftype) bval = self.lspace.basis(bc) #(NC, NQ, ldof//3) @@ -260,7 +260,7 @@ def face_basis(self, bc, index=_S): f2v = self.face_basis_vector(index=index)#(NF, ldof, GD) NF = len(f2v) shape = bc.shape[:-1] - val = bm.zeros((NF,) + shape+(ldof, GD), dtype=self.ftype) + val = bm.zeros((NF,) + shape+(ldof, GD), device=self.device, dtype=self.ftype) bval = self.lspace.basis(bc) #(NF, NQ, ldof//3) @@ -293,7 +293,7 @@ def face_basis_vector(self, index=_S): NF = len(f2e) - f2v = bm.zeros((NF, ldof, GD), dtype=self.ftype) + f2v = bm.zeros((NF, ldof, GD),device=self.device, dtype=self.ftype) # f2v[:, :ldof//2] = e2t[:, 0, None] #(NF, ldof//2, 3) # f2v[:, ldof//2:] = e2n[:, 0, None] f2v = bm.set_at(f2v,(slice(None),slice(None,ldof//2)),e2t[:, 0, None]) @@ -340,7 +340,7 @@ def face_dof_vector(self, index=_S): NF = len(f2e) - f2v = bm.zeros((NF, ldof, GD), dtype=self.ftype) + f2v = bm.zeros((NF, ldof, GD), device=self.device,dtype=self.ftype) # f2v[:, :ldof//2] = e2t[:, 0, None] #(NF, ldof//2, 3) # f2v[:, ldof//2:] = e2n[:, 0, None] f2v = bm.set_at(f2v,(slice(None),slice(None,ldof//2)),e2t[:, 0, None]) @@ -366,7 +366,7 @@ def curl_basis(self, bc): c2v = self.basis_vector()#(NC, ldof, GD) sgval = self.lspace.grad_basis(bc) #(NQ, NC, lldof, GD) - val = bm.zeros((NC,)+(bc.shape[0], )+(ldof, GD), dtype=self.ftype) + val = bm.zeros((NC,)+(bc.shape[0], )+(ldof, GD),device=self.device, dtype=self.ftype) # val[..., :ldof//3, :] = bm.cross(sgval, c2v[:,None, :ldof//3, :]) # val[..., ldof//3:2*(ldof//3), :] = bm.cross(sgval, c2v[:,None, ldof//3:2*(ldof//3), :]) @@ -544,7 +544,7 @@ def set_dirichlet_bc(self, gD, uh, threshold=None, q=None): uh[face2dof[:, ldof//2:]] = bm.sum(gval*f2v[:, ldof//2:], axis=-1) # uh[face2dof[:, :ldof//2]] =0 # uh[face2dof[:, ldof//2:]] =0 - isDDof = bm.zeros(gdof, dtype=bm.bool) + isDDof = bm.zeros(gdof, device=self.device,dtype=bm.bool) isDDof[face2dof] = True return uh,isDDof @@ -564,7 +564,7 @@ def set_neumann_bc(self,gD): points = mesh.bc_to_point(bcs)[isbdFace] n = mesh.face_unit_normal()[isbdFace] hval = gD(points,n) - vec = bm.zeros(gdof, dtype=self.ftype) + vec = bm.zeros(gdof, device=self.device,dtype=self.ftype) k = bm.einsum('fqg, fqlg,q,f->fl', hval, bphi,ws,fm) # (NF, ldof) bm.add.at(vec,face2dof,k) #bm.scatter_add(vec,face2dof,k)